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1.  INTRODUCTION 

Because  today's  industrial  society  consumes  energy  at  such 
an  alarmingly  high  rate,  domestic  supplies  of  fossil  fuels  such  as 
natural  gas  and  oil,  which  are  non-renewable  resources,  have  become 
severely  depleted.   This  high  consumption  rate  is  the  cause  of  the 
current  energy  dilemma  facing  the  world.   The  skyrocketing  cost  of 
these  fossil  fuels  has  forced  a  re-evaluation  of  raw  energy  sources 
used  to  generate  electricity.   Nuclear  fuel  is  a  possible  candidate 
to  replace  fossil  fuel  for  large-scale  power  generator.   However,  public 
fear  of  nuclear  energy  and  its  consequences  has  proved  to  be  a  formidable 
barrier  to  the  development  of  nuclear  power.   Hence,  other  means  to  gener- 
ate power  have  been  investigated.   To  many,  renewable  energy  resources 
such  as  solar  and  wind  energy  are  especially  attractive.   In  addition,  to 
low  input  energy  costs  afforded  by  solar  and  wind  energy  systems,  these 
resources  are  available  to  all  countries  and  consequently  need  not  be 
imported  from  other  lands.   Such  resources  could  help  make  the  U.S. 
and  other  industrialized  nations  less  dependent  on  expensive  foreign 
oil. 

Potentially  available  wind  power,  a  form  of  solar  power,  in 
the  Continental  U.S.,  the  Aleutian  Arc,  and  the  Eastern  Seaboard 
has  been  estimated  to  be  10  GW  [1].   This  is  more  than  thirty  times 
the  estimated  total  annual  average  power  requirements  of  the  U.S.  by 
1980. 

Because  of  this  large  potential  and  the  recent  large  increase 
in  traditional  energy  costs,  the  development  of  various  wind  energy 


systems  is  currently  receiving  vide  attention.   However,  wind  power 
is  not  a  recent  development,  but  has  been  used  by  man  throughout 
history.   As  far  back  as  200  B.C.  wind  turbines  were  used  to  grind 
grain.   From  then  up  until  the  present  day,  wind  turbines  have  been 
used  by  many  societies  to  supply  energy  for  a  variety  of  industrial 
activities.   In  the  mid-19th  century  wind  power  was  used  in  rural 
America  to  pump  water  and  generate  electricity  until  the  1930 's  when 
the  Rural  Electrification  Administration  displaced  these  units  with 
more  reliable  centralized  electric  power.   Thus,  small  wind  systems 
became  a  second  choice  for  power  generation. 

Past  wind  systems  have  not  been  limited  to  serving  small  local 
power  requirements.   A  central  station  wind-powered  electrical  system 
was  conceived  by  Palmer  Putnam  and  built  by  the  S.  Morgan  Smith  Company 
of  York,  Pennsylvania  during  the  early  1940s  [1].   This  system  was 
rated  at  1250  kW  and  fed  power  into  the  electrical  network  of  the 
Central  Vermont  Public  Service  Company.   The  unit  operated  inter- 
mittently until  March,  1945.   At  this  time,  one  of  the  blades  broke 
off  near  the  hub  where  a  known  weakness  had  been  identified  but 
had  not  been  corrected  because  of  the  wartime  material  shortages. 
However,  after  a  comprehensive  economic  analysis,  it  was  discovered 
that  the  unit,  even  if  repaired,  could  not  compete  effectively,  at 
that  time,  with  conventional  electrical  generation  plants.   Hence, 
although  the  project  was  a  technical  success,  economics  forced  the 
plant  to  be  abandoned.   Other  similar  endeavors  in  Denmark,  Russia, 
France,  England,  and  Germany  during  this  same  period  of  time  suffered 
the  same  economic  fate  [1]. 


Currently,  the  rising  cost  of  traditional  energy  sources  has 
kindled  new  interest  in  wind  power.   In  September  1975  under  a  Federal 
wind  energy  program,  the  Energy  Research  and  Development  Administration 
(ERDA)  began  testing  a  wind  turbine  generator  rated  at  100  kW  in  18  mph 
winds  (the  Plum  Brook  unit  of  the  National  Aeronautics  and  Space  Admin- 
istration (NASA)).   Through  this  program,  the  problems  associated  with 
large  wind  turbine  generators  are  to  be  investigated.   These  include  re- 
ducing capital  costs,  eliminating  television  interference,  and  improving 
aesthetic  appearances  [2],   Although  large  wind  turbines  may  offer  the 
potential  of  lower  capital  costs  per  installed  kilowatt,  small  wind  tur- 
bines of  comparatively  simpler  design  are  already  on  the  market  and  do 
not  face  the  problems  of  television  interference  or  aesthetic  impact. 
Consequently,  wind  power  generation  faces  an  economy  of  scale,  i.e., 
whether  to  install  one  large  central  station  unit  or  use  several  small 
units  to  serve  particular  load  demands. 

Whichever  size  units  are  chosen  to  be  used,  areas  of  high  wind 
power  potential  must  be  identified.   Reed  [3]  has  compiled  extensive 
power  calculations  from  historical  meteorological  data  tapes  for  many 
locations.   However,  because  of  the  wind's  inherent  variability,  the 
power  output  from  wind  turbines  is  often  highly  variable  from  day  to 
day  and  even  from  minute  to  minute.   Consequently,  there  have  been  many 
studies  to  try  to  characterize  this  variability  analytically  so  as  to 
be  able  to  predict  the  energy  that  can  be  extracted  by  a  particular  wind 
turbine.   Justus,  et.  al. ,  [4],  made  use  of  the  Weibull  distribution  in 
modeling  wind  speed  distributions  and  computing  total  energy  production 


for  two  different  central  station  power  units,  a  100  kW  wind  turbine 
with  characteristics  similar  to  the  NASA  experimental  unit  and  a  hypo- 
thetical 1  MW  unit,  at  various  sites  throughout  the  United  States. 
Johnson  [5],  too,  used  the  Weibull  distribution  in  analyzing  wind  turbine 
performance  at  specific  Kansas  locations.   In  addition,  Hennessey  [6]  in- 
corporated the  Weibull  distribution  in  his  study  of  computing  mean  power 
densities  at  several  locations;  but  he  solved  for  the  Weibull  distribution's 
parameters  using  a  matching-moments  method,  whereas  Justus  and  Johnson 
estimate  the  parameters  by  a  linear  least  squares  technique.   Corotis  [7] 
also  investigated  the  Weibull  distribution  using  the  matching-moments 
parameter  estimation  method  for  computing  the  available  power  in  the 
wind  as  a  function  of  wind  speed.   He  also  examined  the  use  of  the 
Rayleigh  distribution,  which  is  a  special  case  of  the  Weibull  distri- 
bution involving  the  estimation  of  only  one  parameter.   Cliff  [8] 
also  employed  the  Rayleigh  distribution  to  estimate  the  average  annual 
output  of  a  wind  turbine.   Finally,  Kaminsky  [9]  analyzed  four  different 
distribution  functions,  log-normal,  gamma,  Weibull,  and  Rayleigh,  and 
compared  how  well  each  of  these  distributions  fit  given  wind  speed  data. 
Kaminsky  solved  for  the  parameters  of  the  above  distributions  by  using 
maximum  likelihood  estimators  for  all  except  the  log-normal  distribution 
for  which  he  employed  a  matching-moments  estimation.   However,  common  to 
all  these  studies  that  analyze  wind  turbine  performance  by  using  an  analytical 
distribution  to  characterize  wind  speed  data,  it  is  assumed  that  all  the 
generated  electricity  is  used  by  the  electrical  network,  i.e.,  the  wind 
turbine  is  a  replacement  for  conventional  generators  and  no  generated 
energy  is  wasted. 


Obermeier  [10]  studied  the  prospect  of  using  wind  power  to 
meet  particular  load  demands.   He  investigated  the  use  of  wind 
turbines  in  current  commercial  production  and  his  criterion  for  an 
"optimal"  wind  turbine  and  battery  size  was  one  which  will  supply  the 
entire  demand,  i.e.,  the  customer  does  not  need  to  purchase  any  electri- 
city.  However,  there  was  no  consideration  of  the  cost  to  the  customer 
of  such  an  optimal  system.   Because  economics  plays  a  key  role  in  the 
determination  of  energy  policies,  another  criterion  for  choosing  an 
"optimal"  wind  turbine  size  would  be  to  test  many  different  sized  wind 
systems  and  see  which  one  would  save  the  user  the  most  money  in  com- 
parison to  purchasing  all  the  demanded  electricity  from  a  utility. 

Based  upon  these  previous  wind  power  and  feasibility  studies,  the 
scope  of  this  work  was  two-fold.   First,  various  analytical  models 
of  wind  speed  distributions  were  investigated.   Several  of  the  more 
common  methods  of  solving  for  the  parameters  of  the  Weibull  distribution 
were  studied.   In  addition,  the  beta  distribution  was  introduced  and  in- 
vestigated to  determine  how  well  it  could  represent  wind  speed  data. 
Two  goodness  of  fit  tests  were  performed  on  each  analytical  model  to 
determine  the  appropriateness  of  each  model  in  describing  observed  wind 
speed  distributions.   This  study  of  the  representation  of  observed  wind 
speed  distributions  by  analytically  fit  functions  is  presented  in 
Chapter  2.   In  the  second  phase  of  this  work  a  methodology  was  developed 
whereby  an  appropriate  analytical  wind  speed  model  was  used  to  compute  an 
economically  optimal  wind  turbine  system  to  serve  a  particular  load. 


The  details  of  this  optimization  methodology  and  an  investigation  of 
the  sensitivity  of  the  optimally  sized  wind  turbine  generator  system 
to  the  problem  parameters  is  presented  in  Chapter  3. 


2.   ANALYTICAL  REPRESENTATIONS  OF  WIND  SPEED  DISTRIBUTIONS 

2 .1   Introduction 

To  evaluate  accurately  the  energy  potential  of  a  wind  turbine 
generator  system  (WTGS)  for  a  given  application  and  geographical 
area,  it  is  first  necessary  to  have  sufficient  information  about  the 
wind  speeds  likely  to  be  encountered  by  a  WTGS  at  a  specific  location. 
Wind  speed  characteristics  may  vary  widely  for  different  geographical 
regions  as  well  as  locally  as  a  result  of  local  terrain  features.   More- 
over, the  wind  speeds  vary  throughout  the  day  as  well  as  with  the  season. 
Ideally,  one  would  like  to  have  the  average  wind  speed  distribution  for 
times  throughout  a  day  and  for  days  throughout  the  year  at  a  specific 
site.   For  more  detailed  analyses,  more  detailed  information  about  the 
wind  speed  characteristics  may  be  needed,  e.g.,  standard  deviations  of 
the  average  wind  speed,  wind  duration  (or  persistence)  characteristics, 
and  gust  characteristics.   Such  characterizations  of  the  wind  speeds  may 
be  obtained  from  historical  wind  velocity  data  which  have  been  collected 
at  many  locations  in  the  United  States  and  other  countries  over  a 
period  of  many  years.   In  their  most  basic  form,  these  data  represent 
measured  values  of  wind  speeds  (either  "instantaneous"  or  averaged) 
at  periodic  times  throughout  the  day  for  every  day  of  the  year.   Be- 
cause of  the  large  amount  of  information,  these  historical  wind  records 
are  most  often  stored  on  magnetic  tape  which  can  readily  be  processed 


by  computers.   From  the  historical  records  a  great  deal  of  information 
about  the  average  characteristics  of  the  wind  at  the  specific  location 
can  be  computed  if  a  sufficient  number  of  years  has  been  included  in 
the  data  records. 

For  the  evaluation  of  a  WTGS  at  a  specific  site,  the  most  fund- 
amental wind  speed  information  needed  is  the  average  speed  distribution 
for  various  times  throughout  a  day  in  any  season  (or  month) .   From  the 
meteorological  wind  data  tapes,  various  wind  speed  frequency  (or 
probability)  distributions  have  been  obtained  by  averaging  the  observed 
wind  speed  data  from  various  intervals  throughout  the  day  for  each 
month.   Extensive  compilations  of  such  frequency  distributions  for 
many  locations  in  the  United  States  have  been  published  [11,12],   With 
such  averaged  wind  speed  distributions  the  potential  of  a  given  site 
for  generating  power  can  be  evaluated  accurately. 

However,  even  the  use  of  these  historical  averaged  wind  speed 
distributions  (usually  presented  as  histograms)  still  involves  exten- 
sive data  manipulation.   Moreover,  if  the  average  speed  distributions 
have  been  generated  from  only  a  few  years  of  historical  data,  the 
resultant  frequency  distribution  may  still  differ  significantly  from 
the  actual  wind  frequency  distribution  at  the  site  as  a  consequence 
of  large  wind  speed  variations.   For  the  analysis  of  wind  potential, 
it  is  often  easier  to  use  a  smooth,  analytical,  functional  representation 
of  the  observed  wind  frequency.   The  use  of  such  analytical  functions 
to  represent  the  wind  speed  distribution  allows  a  significant  simpli- 


fixation  in  the  calculation  of  energy  generated  from  a  WTGS  by 
allowing  much  analytical  simplification  in  the  analysis  and,  thereby, 
requiring  considerably  fewer  numerical  computations.   Furthermore,  the 
use  of  analytical  functional  representations  of  the  wind  speed  distri- 
butions tends  to  smooth  out  data  variations  resulting  from  insufficient 
experimental  data.   Strict  use  of  the  experimental  data  may  not  produce 
an  accurate  representation  of  the  actual  wind  speed  distribution  as  a 
result  of  statistical  fluctuations  in  the  data,  whereas  the  analytical 
function  representation  may  yield  a  more  accurate  wind  speed  representation 
by  averaging  over  the  statistical  fluctuations.   Of  course,  the  function 
used  to  represent  the  wind  data  must  be  shown  to  be  very  representative 
of  actual  wind  speed  distributions.   Finally,  most  experimental  wind  speed 
distributions  are  presented  as  histograms  with  relatively  large  velocity 
intervals  or  bins.   This  "binning"  procedure  of  the  historical  measurements 
produces  at  best  an  approximate  representation  of  the  wind  speed  distribution, 
which,  in  reality  is  a  continuous  distribution.   Thus,  the  use  of  analytical 
functions  to  represent  the  wind  speed  distributions  resembles  the  smooth 
and  continuous  nature  of  actual  wind  speed  distributions  as  well  as  pre- 
serving the  distribution  character  and  accuracy  of  representing  the  average 
wind  speeds. 

For  modeling  accurately  power  generation  over  the  speed  range  in 
which  the  WTGS  power  output  varies  rapidly  with  wind  speed,  the  use  of  a 
continuous  wind  speed  model  is  conceptually  more  appealing  than  the  use 
of  a  discrete  frequency  distribution.   The  WTGS  power  output  varies,  ideally, 
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with  the  cube  of  the  wind  speed.   The  discrete  frequency  distribution 
of  wind  speeds  obtained  from  historical  wind  data  usually  contains  only 
two  or  three  intervals  over  the  rapidly  varying  (with  the  cube  of  the 
wind  speed)  transition  range  (from  cut-in  to  rated  speeds)  of  the  WTGS 
response. 

Consequently,  there  has  been  much  effort  in  the  past  few  years 
to  find  various  simple  functions  which  can  represent  accurately  actual 
wind  speed  distributions  as  well  as  to  develop  techniques  for  finding 
the  parameters  of  such  analytical  representations,  which  fit  the  observed 
wind  data  [4,5,6,7,8,9,12].   In  the  first  phase  of  this  work,  several 
techniques  for  fitting  two  functions  (the  Weibull  and  beta  distributions), 
which  are  well  suited  to  represent  wind  speed  distributions,  are  ex- 
amined and  two  tests  (chi-square  and  power  ratio),  are  developed  to 
indicate  the  accuracy  of  the  resulting  fits. 

2 . 2  Description  of  Wind  Speeds 

For  analysis  of  WTGS  performance,  in  addition  to  the  WTTGS  response 
to  wind  speeds,  only  information  about  the  distribution  of  speeds  for 
a  given  daily  time  interval  is  needed.   The  directional  dependence  of 
wind  speeds  is  not  a  concern  in  this  study.   Wind  speed  distributions 
can  be  described  by  either  of  the  following  two  distributions: 

(i)   f(v)  =   "probability  density  function".   The  quantity 
f(v)dv  is  the  probability  that  the  wind  is  in 
the  interval  dv  about  speed  v  and  is  normalized 
such  that 


f(v)dv  =  1.  (2.2-1) 

0 


11 


H 

L)   F(v)  E  J  f(v')dv'  = 


(li)  F(v)  =  J  f(v')dv'  »  "cumulative  distribution"  function, 
i.e.,  the  probability  that  the  wind 
speed  is  less  than  v.  (2.2-2) 


Historical  wind  speed  recordings  are  not  continuous,  but  rather 

the  data  are  grouped  into  discrete  speed  intervals  called  wind  speed 

bins.   In  such  a  discrete  representation,  the  possible  range  of  wind 

speeds,  from  zero  to  high  storm  speeds,  is  divided  into  n  contiguous 

subintervals  of  widths  Av.  bounded  by  speeds  v.  and  v.,,.   In  most 
1  l      l+l 

compilations  of  wind  speed  data,  the  frequency  distribution   of  wind 

speeds  in  each  discrete  subinterval  is  given,  i.e.,  the  probability 

P.  (or  fraction  of  the  total  number  of  observations)  in  which  the  wind 

speed  was  observed  to  be  in  the  i-th  speed  subinterval  (v.<v<v  ,),  is 

i    i+1 

given:  . 


1+1  f(v)  dv.  (2.2-3) 


An  example  of  a  frequency  distribution  compiled  from  historical  records 
is  given  in  Table  2.2-1.   From  such  frequency  data,  a  discretized  form 
of  the  probability  density  function  can  be  constructed  as 

P. 

f .  «  T—  ,  v.  <  v  <  v.,_  .  (2.2-4) 

x   Av.     l        l+l  K     ' 

l 

Also,  the  wind  speed  data  can  be  represented  by  a  discrete  form  of 
the  cumulative  distribution  function, 


*,  =  I     V  J-l,..., n.  (2.2-5) 

i=l 

The  corresponding  probability  density  and  cumulative  distribution 
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functions  of  the  data  in  Table  2.2-1  are  shown  in  Figs.  2.2-1  and 

2.2-2,  respectively.   The  mean  wind  speed,  u,  and  the  dispersion  (or 

2 
variance) ,  o  ,  of  the  wind  speeds  at  any  particular  time  are  calculated 

from  the  continuous  speed  distributions  as 

V   »  I   vf(v)dv,  (2.2-6) 

0 
and 

(v-u)2f(v)dv.  (2.2-7) 


0 
From  the  discrete  representation,  i.e.,  the  experimental  data,  the 

mean  and  variance  can  be  estimated  by 


1=1 


and 


s2  =  I      (v  ...  -  v)2P  (2.2-9) 

1=1   ***       X 

where  n  is  the  total  number  of  speed  subintervals,  and  v...  is  the  mid- 

point  of  wind  speed  subinterval  i. 


2.3  Weibull  Probability  Density  Function 

The  two-parameter  Weibull  distribution  [4,5]  is  the  most  widely 
used  distribution  to  characterize  wind  data  because  of  the  simplicity 
of  its  cumulative  distribution  function.   The  Weibull  distribution  is 
given  as 
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Fig.  2.2-1.   Probability  density  function  for  data  in  Table  2.2-1. 
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Fig.  2.2-2.   Cumulative  distribution  function  for  data  in  Table  2.2-1. 
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f(v)  =   Cf)(f)k_1  exp(-(^)k),    v,k,c>0,  (2.3-1) 

where  v  is  the  wind  speed,  c  is  the  scale  parameter   and  k  is  the 
shape  parameter.      The  corresponding  cumulative  distribution  is 
given  by 

/Vf(v)  dv  =  F(v)  -  1  -  exp[-(-)k].  (2.3-2) 

The  mean  and  variance  of  the  Weibull  distribution,  which  will 
be  needed  later  for  fitting  wind  data,  are  [4]: 

H  -«*(!+£),  (2.3-3) 

and 

2     2      2      2  ?      i 
a     =  c  r(l  +£)  -  cV(l  +i),  (2.3.4) 

where  r(x)  is  the  gamma  function,  defined  by 

r(x)  =  j  y    e  y  dy,   x>0.  (2  3-5) 

0 

The  effect  of  variation  in  the  shape  parameter,  k,  upon  f(v)  is 

illustrated  in  Fig.  2.3-1.   In  this  case  the  scale  parameter,  c,  is 

set  to  unity;  however,  f (v)  can  be  obtained  for  other  values  of  c 

from  these  graphs  by  simply  dividing  the  ordinate  by  c  and  multiplying 

the  abscissa  by  c.   This  adjustment  preserves  the  requirement  of  unit 

area  under  the  curve  [5].   For  values  of- k  between  zero  and  unity,  the 

distribution  has  a  mode  at  zero  and  is  monotonically  decreasing  (ex- 
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2.3-1.    The    standard    form  of    the  Weibull 
distribution. 
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potential-shaped).   For  k  equal  to  unity,  the  distribution  is 

purely  exponential.   For  v  greater  than  3.5,  the  distribution  becomes 

approximately  normal  [6]. 

The  Rayleigh  distribution  is  a  special  case  of  the  Weibull 
distribution.   In  the  Rayleigh  distribution,  the  shape  parameter,  k, 
is  fixed  at  a  value  equal  to  2.   Consequently,  the  Rayleigh  distribution 
is  a  single  parameter  distribution  depending  solely  upon  the  value  of 
the  mean  of  the  distribution. 

2.3.1  Representation  of  Wind  Speeds  by  the  Weibull  Distribution 

(a)     Double  Logarithmic  Transformation  Least  Squares    (DLTLS)  Method 
Several  methods  to  estimate  the  unknown  parameters,  c  and  k,  of  the 
Weibull  distribution  have  been  proposed  [A ,5 ,6, 7 , 12] .   One  method  is  to 
perform  a  least  squares  fit  of  the  observed  data  to  the  double  logarithm 
of  the  observed  cumulative  distribution  function.   The  cumulative  Weibull 
distribution  function  is  a  more  tractable  form  than  its  probability  density 
function.   The  cumulative  distribution  function  is  linearized  by  taking 
the  logarithm  twice  of  each  side  as  shown  below  [5] 

ln[-ln(l-F.)]  =  k  lire.^  -  klnc  .  (2.3-6) 

This  result  is  of  the  linear  form 

y  =  ax  +  b,  (2.3-7) 
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where 

y  -  ln[-ln(l  -  F  )], 

a  =  k, 

and  b  =  -k  lnc 

Hence,  a  least  squares  procedure  can  be  used  on  Eq.  (2.3-7)  to 
yield  estimates  of  the  parameters  a  and  b,  which  in  turn  can  be  used 
to  form  estimates  of  k  and  c.   Thus,  use  of  the  observed  cumulative 
distribution  function,  F.,  with  a  standard  linear  least  squares 
method  gives  estimates  of  the  parameters  c  and  k.   To  avoid  the  case 
where  F.  is  exactly  unity  (i.e.,  when  j  is  equal  to  n) ,  a  modification 
of  Eq.  (2.3-6)  has  been  used  [13],  namely 

ln[-ln(1.0001-F.)]  =  k  lnv.,,  -  k  Inc.  (2.3-8) 

J  i"f*5 

Because  of  the  addition  of  the  small  amount  in  the  left  hand  side 
of  Eq.  (2.3-8),  each  value  is  slightly  in  error.   Alternately,  the 
final  data  point,  F  ,  can  be  omitted  and  Eq.  (2.3-8)  used  with  no 
modification. 

The  least  squares  estimates  a  and  b  (Eq.  (2.3-7))  are  given  by  [14]: 


y  w.  7  w.x.y.  -  J  w.x.  V  w.y. 
.-,  i  .'-,  i  i/i  .L,  i  i  .'.  11 
i=l    i=l 1=1 1=1 

n      n    ~     n      9 

y  w.   y  w.x.  -  (  y  w.x. )" 
■  i   i  .  ,  1 1      .  ,  1 1 

1=1    i=l       i=l 


(2.3-9) 


and 


i=l   l      1-1 


1-1  x  i=l  x 


(2.3-10) 
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where  w.  is  the  weighting  factor  for  speed  subinterval  i. 

Often  the  weighting  factors  for  the  least  squares  technique 
are  assumed  to  be  unity,  i.e.,  each  point  has  the  same  effect  on 
the  least  squares  estimators.  Johnson  [13]  reports  that  weighting 
factors  of  unity  cause  all  speed  subintervals  to  have  the  same  effect 
on  the  least  squares  estimators,  whereas  some  speed  subintervals  actually 
represent  more  observations  than  others.   Consequently,  in  order  for  the 
speed  subintervals  with  the  most  observations  to  have  more  effect  on  the 
least  squares  estimators,  each  speed  subinterval  can  be  weighted  with  the 
probability,  P.,  of  a  wind  speed  occurring  in  that  i-th  speed  sub- 
interval.   Thus,  the  resulting  Weibull  distribution  will  fit  the 
speed  subintervals  which  have  the  most  observations  better  than  the 
Weibull  distribution  obtained  by  using  unit  weighting  factors.   How- 
ever, this  weighting  scheme  is  strictly  judgmental  and  not  statistically 
justified. 

Weighting  factors  which  are  inversely  proportional  to  the  variance 
of  the  residuals  between  the  observed  data  and  the  expected  value  (as 
predicted  by  the  functional  form  used  to  model  the  data)  yield  least 
squares  estimators,  which  are  of  minimum  variance  [15].   Hence,  these 
are  the  best  estimators,  in  the  least  squares  sense.   The  weighting 
factors  are  not  independent  of  the  parameters  of  the  model;  however, 
it  is  normally  assumed  that  the  variance  of  the  residual  is  caused  only 
by  the  variance  in  the  observed  data.   Thus,  the  weighting  factors 
are  commonly  taken  as  the  reciprocal  of  the  variance  of  the  observed 
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data.   In  this  case  (finding  least  squares  estimates  for  a  and  b  of 

Eq .  (2.3-7)),  the  weighting  factors  should  be  the  reciprocals  of  the 

variance  of  the  y  values,  i.e.,  the  reciprocal  of  the  logarithm  of 

the  negative  logarithm  of  unity  minus  the  cumulative  distribution 

of  the  observed  data!   However,  because  the  variance  of  the  observed 

data  is  unknown,  the  statistically  correct  weighting  scheme  is  unable 

to  be  performed. 

Irrespective  of  which  weighting  factors  are  used,  the  above 

linear  least  squares  method  minimizes  the  sum  of  the  errors  associated 

with  the  linear  approximation  of  the  doubly  logarithmic  transformed 

cumulative  distribution  function,  i.e.,  the  fit  parameters  are  chosen 

so  as  to  minimize 

n  9 

E.  =     {(k  lav..,  -  k  lnc)  -  ln[-ln(l  -  F.)]}  .      (2.3-11) 

1      .=1  m  1 

However,  these  methods  do  not  guarantee  that  the  sum  of  the  squared 
errors  of  the  actual  cumulative  distribution, 

EQ  =  I      {(1  -  exp[-[-i^]k])  -  F.}2,  (2.3-12) 

i=l 

is  a  minimum. 

To  illustrate  this  point,  suppose  the  wind  data,  which  are 
shown  in  Table  2.3-1,  are  used.   A  double  logarithmic  transformed  least 
squares  fit  to  the  cumulative  Weibull  distribution  yields  k  equal  to 
1,42  and  c  equal  to  6.20  so  that  the  fitted  distribution  can  be  written  as, 


0-2 
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0.092 
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0.2288 
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0.441 
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21 

0.0624 
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f(v)  =  0.229  (T^o)0'42  expt-(^)1-42]  .  (2.3-13) 


Table  2.3-1.   Hypothetical  Wind  Speed  Data  to  be  Fit  by  a  Weibull 
Distribution  Using  the  Method  of  Least  Squares. 

Speed       Subinterval 

Subinterval    midpoint         P.  f.         F. 

i  i  i 

(knots)        (knots) 


0.0460  0.092 

0.0572  0.3208 

0.0735  0.7618 

0.0293  0.9376 

0.0104  1.0 


The  error  parameter  E. ,  using  this  result  in  Eq.  (2.3-12),  is 

-3 
6.8  x  10   ,  while  E. ,  i.e.,  the  sum  of  the  squared  differences  between 

the  doubly  logarithmic  transform  of  the  cumulative  distribution  associated 
with  Eq.  (2.3-13)  and  the  given  data,  is  0.491.   However,  a  Weibull 
distribution  can  be  found  which  passes  through  all  the  data  points 
except  the  third  speed  subinterval.   The  parameters  of  this  distri- 
bution are  k  equal  to  1.5  and  c  equal  to  10.0  while  written  explicitly 

f(v)  =  0.15  C^)0'5  exp[-  (ygfo)1'5]  •  (2.3-14) 

The  error  parameter,  EL,  using  this  result  in  Eq.  (2.3-12),  is  now 

-3 
smaller  with  a  value  of  0.166  x  10   ,  while  E  increases  to  3.08.   The 

fit  of  Eq.  (2.3-14)  to  the  data  is  considerably  improved,  although  that 
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of  the  doubly  logarithmic  transform  of  the  cumulative  Weibull  is  poorer 
than  in  the  first  case  (Eq.  (2.3-13)).   Hence,  the  DLTLS  estimation  de- 
scribed in  this  section  assures  the  minimization  of  E  ,  but  not 
necessarily  the  minimization  of  Eq ! 

(b)     Hatching-Moments  Method 

Another  procedure  to  obtain  estimates  of  c  and  k  is  by  a  matching- 
moments  method.   The  sample  mean  and  variance  of  the  wind  data,  given  by 
Eqs.  (2.2-8)  and  (2.2-9),  are  set  equal  to  the  mean  and  variance  of  the 
Weibull  probability  density  function,  Eqs.  (2.3-3)  and  (2.3-4),  to 
obtain 

v  =  V      =  cr(l  +  h,  (2.3-15) 

k 

s2  =  a2   =  c2[F(l  +  |)  -  r2(l  +  i)J.  (2.3-16) 

From  these  equations,  c  can  be  expressed  in  terms  of  k  as 

(2.3-17) 


r(i  +  i) 


Substitution  of  c  from  Eq.  (2.3-17)  into  Eq.  (2.3-16)  yields  an  equation 
for  k,  namely, 

r  a  +  h      -2      2 

v  I  s  =  0.  (2.3-18) 


r2(i  +  i) 


To  solve  Eq.  (2.3-18)  for  k,  a  numerical  procedure  for  finding 
the  roots  of  a  nonlinear  equation  must  be  used.  In  this  study,  the 
Mueller's  iteration  method  (an  elegant  successive  bisection  technique  [30]) 


24 


is  preferred  over  a  technique  such  as  Newton's  method,  because  the 
former  does  not  require  the  calculation  of  derivatives  which  can  become 
cumbersome  when  dealing  with  gamma  functions.   Once  k  is  found,  c  can 
be  determined  by  substitution  of  k  into  Eq.  (2.3-17).   If  calculations 
are  to  be  performed  by  hand,  a  crude  estimate  of  k  can  be  made  by  using 
Kotel'nikov's  nomogram  [6].   The  value  for  c  can  then  be  obtained  by 
substitution  of  k  into  Eq.  (2.3-17). 

The  matching-moments  technique  is  very  flexible  in  that  only 
the  sample  data  mean  and  variance  are  needed  to  solve  for  k  and  c. 
Consequently,  diurnal  effects  on  the  wind  distribution  can  be  studied 
with  some  ease.   Since  diurnal  variations  are  generally  given  in  terms 
of  their  effect  on  mean  seasonal  wind  speeds,  the  resulting  frequency 
distribution  is  easily  computed.   However,  a  numerical  method  or  Kotel'- 
nikov's nomogram  approximation  is  needed  to  solve  for  k;  this  numerical 
computation  is  an  unattractive  feature  of  this  method. 


2.4   Beta  Probability  Density  Function 

The  second  probability  density  function  proposed  for  the  modeling 


}f  wind  data  is  the  two-parameter  beta  distribution, 


f(v) 


e^*-1  (i 


v    B(ct,6)   v 
max         max 

0,    otherwise 


v  .e-i 


0<v<v 
max 


(2.4-1) 


where 


Ka,B)  = 


r(a)r(B) 
:(o  +  B) 
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a  and  &  are  positive  parameters,  and  v    is  the  maximum  speed  for 
r  max  r 

which  the  beta  distribution  is  defined,  i.e.,  maximum  observed  wind 

speed. 

Variation  of  the  parameters  a  and  6  causes  the  beta  distribution  to 

assume  many  shapes  as  illustrated  in  Fig.  2.4-1.   For  a  less  than  £, 

the  beta  distribution  is  positively  skewed,  i.e.,  there  is  a  longer 

"tail"  to  the  right  of  the  maximum  than  to  the  left.   For  a  greater 

than  B  the  distribution  is  negatively  skewed  and  for  a  equal  to  6,  the 

distribution  is  svmmetric  about  the  value  ^T   -In  the  special  case  of 

max 

both  a  and  3  equal  to  unity,  the  beta  distribution  is  uniform  over 
its  entire  defined  wind  speed  range. 

The  mean  and  variance  of  the  beta  distribution,  which  will  be 
needed  in  the  following  subsection  for  fitting  wind  data,  are 

y  =  —TS    v  -   .  (2.4-2) 


(2.4-3) 


(a  +  6)2(a  +  B  +  1)   ™aX 


2.4.1  Representation  of  Wind  Speeds  by  the  Beta  Distribution 

To  use  the  beta  distribution  to  model  given  wind  speed  data,  a 
method  must  be  found  to  obtain  values  for  the  parameters  a  and  3  of 
the  beta  distribution.   Techniques  such  as  a  non-linear  least  squares 
or  maximum  likelihood  estimates  can  be  used  to  find  estimates  of  a 
and  3.   However,  the  above  techniques  require  numerical  procedures  in  the 
computation  of  a   and  &. 
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The  standard  form  of  the  beta  distribution. 
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Another  method  to  solve  for  a  and  g  is  by  the  raatching-moments 

technique.   As  was  done  when  solving  for  the  parameters  of  the 

Weibull  distribution  by  matching-moments,  the  sample  data  mean, 

Eq.  (2.2-8),  and  variance,  Eq.  (2.2-9),  are  set  equal  to  the  beta 

distribution's  mean  (Eq.  (2.4-2))  and  variance  (Eq.  (2.4-3)).   The 

resulting  equations  can  be  solved  simultaneously  for  a  and  3  to  yield 

—     v(v    -  v) 
a=_v_  j ma| _  1](  (24_4) 

max       s 
=  !max V_a    _  (2.4-5) 


It  should  be  emphasized  that  a  and  B  can  be  obtained  analytically. 

The  beta  distribution  is  a  particularly  attractive  distribution 
for  wind  modeling  because  it  can  assume  many  shapes  and  because  it  is 
defined  only  over  a  finite  speed  interval.   These  two  properties  help 
the  beta  distribution  to  resemble  the  given  wind  data  more  closely 
since  all  real  wind  frequency  distributions  are  zero  beyond  some  max- 
imum wind  speed,  v   .   Furthermore,  the  parameters  of  the  beta  distri- 
r      max 

bution  are  more  easily  obtained  by  matching-moments  than  those  for  the 
Weibull  fit  since  a  numerical  (iterative)  solution  is  not  required  for 
the  beta  distribution. 

2.5  Description  of  Goodness  of  Fit  Tests 

Once  the  parameters  of  the  desired  distribution  are  computed, 
the  accuracy  of  these  distributions  to  represent  the  wind  speed  data 
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2 
must  be  verified.   Two  tests,  a  chi-squared  (x  )  test  and  a  power 

ratio   test,  to  judge  the  goodness  of  fit  are  developed  in  this  section. 

(a)     y2  Test 

2 
The  first  test  is  the  x  test  which  is  based  on  the  following 

statistic  [16]: 

,    n   [NP.°bE  -  NP.]2 

X2  -  I  -Ag- J—  ,  (2.5-D 

1=1       i 

where  P  °   is  the  observed  frequency  of  occurrence  of  wind  speeds  in 

i 

subinterval  i  obtained  from  recorded  data,  while  P.  is  the  expected 

frequency  as  predicted  by  the  analytical  model  (i.e.,  as  computed 

from  Eq.  (2.2-3).   The  symbol  8  represents  the  total  number  of  all  wind 

speed  observations  used  to  obtain  the  observed  frequency  distribution. 

2 
To  determine  the  significance  of  the  x   test,  the  number  of 

degrees  of  freedom,    equal  to  the  number  of  speed  subintervals, 

n,  minus  the  number  of  different  independent  linear  restrictions 

imposed  on  the  observations,  must  be  determined.   For  the  present  case, 

one  restriction,  i.e.,  the  loss  of  one  degree  of  freedom,  is  due  to  the 

fact  that  the  probability  in  the  last  speed  subinterval  is  determined 

after  the  probabilities  in  the  first  n-1  speed  subintervals  are  known. 

Furthermore,  an  additional  constraint  (loss  of  one  degree  of  freedom) 

results  from  each  independent  parameter,  e.g.,  v  and  s,  which  allows 

the  determination  of  a  and  B,  of  a  distribution  estimated  from  the  data 

[16].   Thus,  the  number  of  degrees  of  freedom,  v,  is  given  by, 

v  =  n  -  1  -  e,  (2.5-2) 
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where  n  is  the  number  of  speed  subintervals  and  e  is  the  number  of 
parameters  estimated  from  the  data. 

The  number  of  degrees  of  freedom  and  the  calculated  value  of 

2 
X  are  used  to  determine  the  level  of  significance,  i.e.,  the 

2 
probability  that  chance  would  allow  a  value  of  x  as  large  or  larger 

than  the  one  calculated.   This  significance  level  can  be  found  from 

either  standard  tables  or  numerically  integrating  the  incomplete  gamma 

function,  which  the  x  distribution  follows.   Thus,  if  the  calculated 

2  2 

X  is  larger  than  the  theoretical  x  >  the  hypothesis  that  the  analytical 

distribution,  with  parameters  estimated  from  the  data,  represents  the 
wind  speed  data  is  rejected. 

In  the  computation  of  x  ,  it  is  imperative  that  NP .     is  the 
number  of  observations  and  not  a  probability  percentage.    The  x 
test  is  invalid  if  probability  percentages  are  used  unless  the  sample 

size  is  exactly  100  [16],   The  sample  size  is  crucial  because  as  noted 

2 
from  Eq.  (2.5-1),  the  x  value  is  a  function  of  the  sample  size.   Hence, 

2 
use  of  the  incorrect  sample  size  changes  the  x  value  and  consequently 

2 
changes  the  significance  level  of  the  x   test  which  could  result  in 

an  erroneous  conclusion  about  the  fit  of"  data  to  the  analytical 

distributions. 

2 
Also,  when  performing  the  x   test,  wind  speed  subintervals  should 

be  grouped  such  that  each  wind  speed  subinterval  for  the  analytical 

distribution  contains  at  least  a  single  observation  [17].   With  this 
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2 
adjustment,  large  contributions  to  x   from  wind  speed  subintervals 

with  few  observations  are  avoided.   Table  2.5-1  illustrates  the 

2 
effect  of  grouping  data.   This  example  shows  that  the  statistic  x 

can  be  reduced  by  a  factor  of  2  if  the  data  are  grouped  properly,  i.e., 

combining  wind  speed  subintervals  so  that  every  subinterval  has  at  least 

a  single  observation.   It  must  be  noted,  though,  that  grouping  the 

intervals  at  the  tails  so  that  the  expected  observations  are  much  greater 

2 
than  unity  causes  the  x~  test  to  loose  its  power.   This  is  because 

grouping  may  cover  up  the  most  distinct  differences  where  the  two 

distributions  differ  the  most  [17]. 

(b)     Power  Ratio  Test 

In  the  second  test,  the  power  available  in  the  wind  for  a  given 
analytical  distribution  is  computed  and  compared  to  the  power  avail- 
able as  computed  from  the  histogram  of  data  from  which  the  analytical 
model  was  obtained.   The  ratio  of  these  two  power  calculations  is  de- 
fined to  be  the  power  ratio.      A  power  ratio  of  unity  indicates 
power  calculations  using  the  observed  wind  speed  distribution  and 
the  fitted  analytical  distribution  yield  identical  results.   Conse- 
quently, this  test  shows  the  accuracy  of  the  analytical  fit  over 
the  entire  speed  range  of  interest  for  wind  turbines  rather  than 
comparing  accuracies  at  selected  intervals. 

To  compute  this  power  ratio,  an  expression  for  average  wind 
power  must  be  obtained.   Using  the  analytical  distribution,  the 
average  power  output  of  a  WTGS,  P.,  ,  is  given  by  the  equation 
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2 
Table  2.5-1.   Tabulation  of  y.  Values  Illustrating  Effect  of  Grouping 

Data. 


Wind  Speed 
Subinterval 

NP0bs 

NP. 

Contribution 
to*2 

Contribution 
to  x2 

(knots) 

(ungroupe 

d) 

(grouped) 

0.0 

-  1.0 

61 

56.3 

0.392 

0.392 

1.0 

-  3.5 

293 

305.8 

0.552 

0.552 

3.5 

-  6.5 

A52 

434.5 

0.747 

0.747 

6.5 

-10.5 

517 

498.6 

0.679 

0.679 

10.5 

-16.5 

452 

455.5 

0.020 

0.020 

16.5 

-21.5 

143 

162.7 

2.45 

2.45 

21.5 

-27.5 

52 

58.8 

1.02 

1.02 

27.5 

-33.5 

10 

6.53 

0.253 

*1.16 

33.5 

-40.5 

2 

0.27 

11.1 

Totals       1982  1982       17.2  7.02 


*last  two  subintervals  combined 
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V         =  j     R(v)f(v)dv,  (2.5-3) 

tit   Q 

where  R(v)  is  the  response  function   of  the  wind  turbine  generator, 
i.e.,  the  power  obtained  from  a  WTGS  when  the  wind  has  a  speed  v, 
and  f(v)  is  the  wind  speed  probability  density  function  to  be  tested. 
To  calculate  the  power  available  when  using  the  discrete  probability 
density  function,  f.,  the  integral  in  Eq.  (2.5-3)  must  be  divided  up 
and  summed  over  each  of  n  speed  subintervals  with  f (v)  equal  to  f.  for 

v  between  v  and  v.  , .   Hence,  when  using  the  discrete  probability 

l      i+l 

density  function  to  compute  the  average  power  output  of  a  WTGS,  p0^s> 
the  following  equation  is  used: 


P  .   =   )  f.  /  i+1  R(v)  dv.  (2.5-4) 

obs    ,**,  i  ' 
1=1   v. 


To  evaluate  the  above  integrals,  the  WTGS  response  function,  R(v), 
must  first  be  selected.   Several  response  function  models  have  been 
proposed.   Justus,  et .  al.,  [4]  suggest  a  quadratic  polynomial.  The 
model  used  for  such  a  fit  was  NASA's  Plum  Brook  WTGS.   The  parameters 
of  this  unit  are  a  rated  power  output  of  100  kW,  a  cut-in  speed   (a 
wind  below  which  the  generator  produces  no  power)  of  8  mph,  a  rated 
speed   of  18  mph,  and  a  furling  speed   (a  wind  speed  above  which  the 
generator  turns  off  to  prevent  damage  to  the  system)  of  60  mph.   The 
response  function  is  described  by 
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0            v<v 
—  c 

R(v)  = 

A  +  Bv  +  Cv",  v  <v<v                   (2.5-5) 
c  —  r 

P  ,              V  <v<v„   , 

r           r  —  furl 

0,           v>v,   , 
furl 

where  v  is  the  cut-in  speed,  v  is  the  rated  speed,  v,  ,  is  the 
c                     r                    furl 

furling  speed,  and  P  is  the  rated  speed.   The  coefficients  A,  B,  and 

C  are  determined  from  the  following  conditions 

2 
A  +  Bv  +  Cv   =0 
c     c 

A  +  Bv  +  Cv   -  P 
r     r     r 

(2.5-6; 

,       v   3 

A  +  Bv  +  Cv   =  P   ( — ) 
o     o     r  v 
r 

J 

where 

v  +  v 
c    r 

V0  "    2     ' 

The  response  function  chosen  for  use  in  this  study  is  a  more 

idealized  one  and  is  based  upon  the  theoretical  power  in  the  wind  due 

to  the  mass  and  velocity  of  air  molecules.   The  total  power  available 

from  the  motion  of  air  with  speed  v  through  a  cross-sectional  area 

A  is  given  by  [5] 

P  =  i  pAv3,                              (2.5-7) 
w   2 

3 
where  P   is  the  power  in  watts,  p  is  the  density  of  air  in  kg/m  , 

w 
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A  is  the  exposed  area  in  m  ,  and  v  is  the  wind  speed  in  m/s.   The 

important  feature  of  this  relation  is  that  the  available  power  varies 

with  the  cube  of  the  wind  speed. 

The  wind  turbine  response  function  used  in  this  study  was  assumed 

to  follow  a  cubic  relation  also.   Fig.  2.5-1  shows  this  idealized  WIGS 

response  function.   Below  v   (cut-in  speed),  the  wind  turbine  produces 

effectively  zero  power  due  to  electrical  and  mechanical  losses. 

Between  v  and  v   (rated  speed),  power  proportional  to  the  cube  of  the 

wind  speed  is  produced.   Once  the  rated  speed  is  reached,  the  wind  turbine 

reaches  its  rated  power  output,  P  .   For  speeds  greater  than  v  ,  a  control 

system  varies  the  blade  pitch  so  that  the  generator  capacity  is  not  exceeded 

and  the  power  ouput  is  maintained  at  its  rated  power.   However,  if  the  wind 

ever  exceeds  v.   , ,  the  WTGS  is  shut  down  or  furled  in  order  to  prevent 
furl  r 

damage  to  the  system.   This  idealized  model  is  expressed  by  the  following 
equation: 


R(v) 


fv  V 


(2.5-8) 


V  <  V  <  V-   -.  , 

r    —  furl 


,  v  >  v 


furl 


Because  storm  level  winds  are  usually  suppressed  in  reported  wind 

speed  data,  the  speed  v,   ,  is  set  equal  to  the  v    which  is  the 
furl  max 

maximum  speed  reported  in  the  wind  speed  data.   Furthermore,  unless 
specified,  the  cut-in  speed  is  given  by  the  equation 
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Fig.  2.5-1. 


Wind  Speed 
An  idealized  KTGS  response  function. 
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v  -  0.46416  v  .  (2.5-9) 

c  r 

This  relation  implies  that  the  WTGS  is  idle  when  the  output  power  is  less 
than  10%  of  the  rated  power  for  the  WTGS. 

Substitution  of  the  WTGS  response  function  given  in  Eq.  (2.5-8) 
into  Eq.  (2.5-3)  yields  an  explicit  form  for  the  average  power  out- 
put of  a  WTGS  using  the  analytical  distribution  as  follows 


Y 

t    =  |  p  rM3f(v 


where  Y  »  min    (v    ,    v        ) . 

1  r        max 

If,  in  any  of  the  integrals  in  Eq.  (2.5-10),  the  upper  limit  is  less 
than  the  lower  limit,  the  integral  is  assumed  to  be  zero.   Because  the 
above  integrals  cannot  be  evaluated  analytically  when  the  analytical 
distribution  is  used,  a  numerical  technique  is  required.   Standard  Gauss- 
Legendre  quadrature  was  used  (see  Section  3.2-4). 

When  using  the  WTGS  response  function  (Eq.  (2.5-8))  to  compute  the 
average  power  output  of  a  WTGS  for  the  discrete  wind  speed  distribution 
from  Eq.  (2.5-4),  care  must  be  taken  to  use  the  proper  form  of  the  WTGS 
response  function  depending  upon  the  magnitude  of  the  speeds  in  a  partic- 
ular speed  subinterval.   For  example,  if  the  rated  speed  occurs  in  the 
middle  of  a  speed  subinterval,  the  speeds  in  the  subinterval  less  than 
or  equal  to  v  follow  the  cubic  power  relationship  given  in  Eq.  (2.5-8) 
while  the  speeds  greater  than  v  follow  the  constant  power  relationship. 
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Consequently,  computation  of  the  power  ratio  requires  calculation 
of  the  average  power  using  both  the  analytical  and  discrete  probability 
density  functions.  The  ratio  of  these  two  results  gives  an  indication 
of  the  closeness  of  fit  between  the  analytical  and  discrete  wind  speed 
probability  density  functions  for  wind  power  calculation.  A  value  of 
unity  means  there  is  an  excellent  fit,  while  values  significantly  less 
than  or  greater  than  unity  indicate  an  under-  or  overestimation  of  the 
available  wind  power,  respectively. 

2.6  Results 

To  solve  for  the  parameters  of  the  fitting  distributions  and  to 
evaluate  the  goodness  of  fit  tests,  a  computer  routine,  CURVEFIT,  was 
developed.   A  listing  and  explanation  are  contained  in  Appendix  A. 
Wind  speed  data  were  obtained  from  the  National  Climatic  Center,  Asheville, 
North  Carolina  for  seventeen  sites  throughout  the  United  States.   A 
map  showing  the  locations  of  the  sites  is  given  in  Fig.  2.6-1.   These 
sites  are  representative  of  the  many  possible  wind  conditions  and  power 
densities  found  throughout  the  United  States.   Table  2.2-1  is  a  typical 
listing  of  the  "binned"  wind  speed  data  obtained  from  the  National 
Climatic  Center.   For  each  of  the  seventeen  sites,  wind  speed  data 
were  given  in  eight,  three-hour  intervals  throughout  a  day  for  a  particular 
month.   Data  for  the  months  of  January,  April,  July,  and  October  were 
chosen  in  this  study  to  simulate  the  four  seasons  of  the  year.   Con- 
sequently, 544  (=  8  x  4  x  17)  observed  wind  speed  distributions  were 
chosen  and  fit  by  analytical  distributions.   Although  anemometer  heights 
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for  the  recorded  wind  speed  data  typically  varied  during  the  recording 
period  at  individual  observation  stations,  no  attempt  could  be  made  to 
adjust  for  these  variations.   This  is  because  complete  historical  records 
of  anemometer  heights  at  each  location  used  were  unavailable.   However, 
anemometer  heights  were  roughly  from  7  to  10  meters  and  consequently, 
any  changes  incurred  because  of  varying  anemometer  heights  should  be 
minimal . 

Example  results  from  the  CURVE FIT  program  are  given  in  Tables  2.6-1 
and  2.6-2.   Two  tables  are  used  to  present  the  results  from  each  location. 
The  first  table  lists  the  parameters  for  each  of  the  analytical  distri- 
butions calculated  from  the  wind  speed  data  while  the  second  table  lists 
the  results  of  the  goodness  of  fit  tests  run  on  each  analytical  distri- 
bution.  In  Table  2.6-1,  the  first  line  gives  the  mean,  standard  deviation, 
and  the  fit  parameters  for  the  wind  speed  distribution  data  from  Syracuse, 
New  York,  in  the  month  of  January  during  the  first  three  hours  of  a 
typical  day  (beginning  at  midnight). 

The  parameters  given  in  Table  2.6-1  can  be  used  to  plot  the 
probability  density  functions.   Figure  2.6-2  shows  the  four  analytical 
fits  whose  parameters  are  given  in  the  first  line  of  Table  2.6-1  to- 
gether with  the  actually  observed  distribution. 

Table  2.6-2  presents  the  goodness  of  fit  statistics  for  the  cases 

2 
presented  in  Table  2.6-1.   This  table  lists  the  x  values  for  each  of 

the  four  analytical  distributions.   The  integer  in  parentheses  is  the 

2 
number  of  degrees  of  freedom  associated  with  the  x  value  for  that 
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Fig.  2.6-2.   Comparison  of  analytical  distributions  to 

observed  distribution  for  Syracuse,  New  York, 
data  (January,  hours  1-3). 
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particular  analytical  fit.   The  final  four  columns  show  the  power 
ratio  computed  for  each  analytical  fit.   The  rated  power  (100  kW) , 
cut-in  speed  (8  mph) ,  and  rated  speed  (18  mph)  used  were  those  of 
the  NASA  experimental  wind  turbine  generator  at  Plum  Brook,  Ohio. 

2 
(a)     Results  of  x  Test 

2 
The  results  of  the  >;  test  reveal  rather  large  values  for 

all  analytical  fits.   For  the  Weibull  distributions  using  the 

2 
least  squares  technique  to  estimate  parameters,  all  X  values 

2 
are  in  excess  of  180.   This  is  quite  large  considering  the  x  value 

at  the  0.995  confidence  level  with  ten  degrees  of  freedom  is  25.2 

2 
and  all  tabulated  x  values  have  fewer  than  ten  degrees  of  freedom. 

Hence,  there  is  less  than  a  0.5%  chance  that  the  least  squares  fitted 

distributions  describe  the  given  data.   This  result  was  not  totally 

unexpected  because,  as  shown  earlier,  the  least  squares  estimation  of 

the  Weibull  distribution  parameters  does  not  guarantee  the  best  fit  to 

the  data. 

2 
However,  large  x  values  were  also  obtained  for  the  matching- 

moments  estimation  of  both  the  Weibull  and  beta  distribution  para- 
meters.  Over  93%  did  not  pass  the  non-significance  test  at  the  0.005 

2 
significance  level.   But  there  are  some  x  values  which  indicate 

a  good  fit  of  the  analytical  functions  to  the  data  (see  Table  2.6-2 

in  which  the  analysis  of  the  fits  to  the  wind  data  from  Syracuse, 

New  York  is  shown).   In  the  examples  of  Table  2.6-2,  about  half  of  the 

2 
match in g-moments  Weibull  distribution  pass  the  x  test  at  the  0.005 

significance  level,  and  of  those  all  except  one  pass  the  test  at  signifi- 
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cance  levels  less  than  or  equal  to  0.01.   For  the  beta  distribution,  only 
31%  pass  the  test  at  a  significance  level  of  0,005  or  less.   Examination 

of  the  remaining  tables  in  Appendix  B  shows  other  isolated  cases  where 

2 

the  x  values  for  both  the  Weibull  and  beta  distribution,  using  the 

matching-moments  parameter  estimation,  yield  significant  results.   In 
addition,  there  are  sporadic  cases  where  the  beta  distribution  gives  the 
largest  x  value  of  all  of  the  analytical  fits. 

Similar  results  were  obtained  by  Kaminsky  for  all  the  analytical 

distributions  he  used  to  describe  observed  wind  speed  data  [9],  i.e., 

2 

the  x  values  were  typically  very  large.   However,  of  the  four  distri- 
butions (log-normal,  gamma,  Weibull,  and  Rayleigh)  he  used  to  characterize 

2 

wind  speed  distributions,  the  lowest  x  values  were  obtained  with  the 

gamma  and  Weibull  distributions. 

2 
Because  the  x   test  merely  sums  the  square  of  the  deviations 

between  the  data  and  the  analytical  fit  at  each  speed  sub interval,  the 

.  2 
contributions  to  the  x  value  from  each  of  these  speed  subintervals  was 

investigated.   In  one  case,  for  the  April  wind  speed  data  from  Columbus, 

Indiana  during  hours  nine  through  twelve,  both  the  Weibull  and  beta 

2 
distributions  using  matching-moments  estimation  yield  x  values  greater 

than  the  other  Weibull  distributions.   A  breakdown  of  the  contributions 

9 

to  the  x  value  by  speed  subinterval  is  shown  in  Table  2.6-3. 

For  the  Weibull  and  beta  distributions  using  matching-moments 

2 
estimation,  98%  of  the  total  contribution  to  the  final  x  value 

arises  from  the  first  two  speed  subintervals,  or  speeds  between  0  and  3.5 
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knots.   The  remaining  speed  subintervals  contribute  little  to  the  final 

X  value.   On  the  other  hand,  the  other  two  Weibull  distributions  (obtained 

by  least  squares  fits)   show  rather  uniform  contributions  from  each  speed 

2 
subinterval  to  the  total  x  value.   There  are  no  speed  subintervals  which 

contribute  such  a  large  amount  to  the  final  value.   This  is  typical  of 

most  of  the  remaining  data,  i.e.,  the  deviations  between  the  Weibull  and 

beta  distribution  using  the  matching-moments  estimation  techniques  and  the 

actual  wind  speed  distribution  occur  in  the  first  two  speed  subintervals, 

whereas  in  the  two  least  square  fits  to  the  Weibull  distribution,  an 

2 
equally  large  amount  is  contributed  to  the  x   statistic  by  all 

speed  subintervals.   Consequently,  this  indicates  that  the  matching- 
moments  estimations  of  the  parameters  of  the  Weibull  and  beta  distri- 
butions fit  the  intermediate  and  high  wind  speed  subintervals  much  more 
closely  than  at  the  low  speed  end  of  the  distribution.   It  is  the  middle 
and  upper  speed  subintervals  which  are  important  in  the  analysis  of  a 
WTGS  and  for  which  the  matching-moments  technique  produces  good  fits 
with  the  Weibull  and  beta  distributions. 

(b)     Results  of  Power  Ratio  Test 

The  power  ratio  values  were  also  obtained  for  all  544  wind  distri- 
butions fit  by  the  Weibull  and  beta  distributions.   In  addition  to  use 
of  the  cut-in  and  rated  speeds  of  the  NASA  Plum  Brook  WTGS  (a  cut-in 
speed  of  8  mph  and  a  rated  speed  of  18  mph) ,  rated  speeds  of  12,  15,  21, 
and  24  mph  were  also  used.   Rated  power  was  held  constant  at  100  kW.   Cut- 
in  speeds  for  the  latter  four  cases  were  given  by  Eq.  (2.5-9). 
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The  mean  power  ratios  obtained  from  the  wind  speed  data  at  the 
seventeen  sites  for  each  WTGS  size  are  tabulated  in  Table  2.6-4.   As 
explained  earlier,  544  observed  wind  speed  distributions  were  used 
to  calculate  the  power  ratios  for  each  WTGS  size.   From  this  table 
it  can  be  seen  that  parameter  estimations  obtained  by  unweighted  and 
weighted  least  squares  of  the  doubly  logarithmic  transformed  cumulative 
Weibull  distribution  grossly  underestimate  the  power  available, 
whereas  the  match in g-moments  parameter  estimations  of  both  the  Weibull 
and  beta  distributions  yield  values  very  close  to  unity. 

Frequency  plots  showing  the  distributions  of  the  power  ratio 
values  about  their  calculated  mean  are  shown  in  Figs.  2.6-3  to  2.6-12. 
Because  of  the  large  dispersion  of  power  ratios  about  the  mean  value 
for  the  Weibull  distributions  which  were  obtained  by  the  least 
squares  estimators,  their  power  ratios  are  grouped  into  wider  class 
intervals  in  order  to  yield  a  smoother  distribution. 

Since  the  power  ratios  obtained  from  both  the  beta  and  Weibull 
distributions  using  matching-moments  parameter  estimation  are  so  close 
to  the  expected  value  of  unity,  it  is  reasonable  to  test  these  values 
and  ascertain  whether  there  is  a  significant  difference.   Consequently, 
a  t  test  [18]  is  applied. 

A  t  test  is  used  to  determine  if  there  is  a  significant  difference 
between  an  observed  mean  and  a  theoretical  mean  or  between  two  observed 
means.   Philosophically,  one  establishes  the  hypothesis  that  the  two  means 
(observed  and  theoretical  or  two  observed  means)  are  equal  or  that  the 


Table  2.6-4-  List  of  Power  Ratio  Means  and  Standard  Deviations  for 
Various  Analytical  Distributions  and  Generator  Rated 
Speeds. 
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Method  For 

Calculation  of 

Power  Ratio 


Rated  Speed  Power  Ratio       Sample 
v  Mean         Std.  Dev. 

(mph)  7  Sr 


Weibull-      1st.      sq. 


Neibull-wt.      1st.      sq. 


Weibull-  matching-inoments 


Beta 


12 

0.589 

0.115 

15 

0.556 

0.109 

18 

0.529 

0.108 

21 

0.519 

0.113 

24 

0.519 

0.118 

12 

0.664 

0.0754 

15 

0.648 

0.0607 

18 

0.636 

0.0730 

21 

0.648 

0.105 

24 

0.671 

0.159 

12 

0.982 

0.0355 

15 

0.988 

0.0376 

18 

0.975 

0.0427 

21 

0.986 

0.0475 

24 

0.998 

0.0647 

12 

0.992 

0.0232 

15 

1.01 

0.0311 

18 

1.001 

0.0331 

21 

1.01 

0.0342 

24 

1.02 

0.047 
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samples  used  to  calculate  the  observed  means  are  from  the  same  pop- 
ulation.  The  calculated  t  value  is  compared  to  a  t  value  (called  the 
critical  t  value)  at  a  specified  level  of  significance,  i.e.,  to  the 
integral  of  the  t  distribution  with  an  integration  limit  set  such 
that  the  area  of  integration  is  equal  to  unity  minus  the  level  of 
significance.   Thus,  the  level  of  significance  is  the  probability  that 
chance  will  allow  a  t  value  equal  to  or  greater  than  the  critical  t 
value. 

The  t  statistic  is  defined  by  [18] 


[sWS 


where  r  is  the  sample  mean,  y  is  the  mean  which  r  is  to  be  tested 
against  (equal  to  unity  for  the  sample  variance,  and  N  is  the  sample 
size),  and  s   is  the  standard  deviation  of  the  sample  values  from  the 
sample  mean. 

Table  2.6-5  lists  the  results  of  this  test  for  the  power  ratio 
sample  means.   The  t  value  at  the  0.001  significance  level  with  an 
infinite  number  of  degrees  of  freedom  (each  mean  power  ratio  mean 
was  computed  from  544  wind  speed  distributions)  is  3.291,  which  means 
that  there  is  less  than  a  0.1%  probability  that  a  random  sample  of 
power  ratios  drawn  from  a  population  with  a  mean  of  u  will  yield  a 
t  value  of  3.291  or  greater.   All  but  two  power  ratio  sample  means 
have  t  values  larger  than  3.291;  thus,  it  can  be  stated  that  the  true 
mean  of  the  power  ratio  using  a  particular  analytical  distribution  is 
something  other  than  unity  at  the  0.001  level  of  significance. 
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Although  the  t  test  does  not  confirm  the  suspicion  that  the  power 
ratio  values  for  the  distributions  using  matching-moments  estimators  have 
a  mean  of  unity,  a  significant  difference  can  be  seen  between  the  power 
ratios  for  the  distributions  using  matching-moments  estimators  and  those 
distributions  using  least  squares  estimators.   The  latter  distributions 
yielded  power  ratio  values  ranging  from  .29  to  2.2  with  most  values 
being  between  .4  and  .8  while  the  former  distributions  were  consistent- 
ly in  the  .85  to  1.4  range  with  90%  of  the  values  between  .90  and  1.1. 
Since  the  speeds  at  which  the  wind  turbine  generators  produce  significant 
power  correspond  to  the  intermediate  and  high  speed  range,  and  since  the 
best  power  ratio  values  are  obtained  from  the  matching-moments  estimation 
of  the  Weibull  and  beta  distribution  parameters,  the  conclusion  drawn  from 
the  x  test  that  the  matching-moments  estimators  yield  distributions  which 
accurately  fit  the  intermediate  and  high  speed  subintervals  is  substan- 
tiated.  Furthermore,  because  a  wind  turbine  cannot  use  low  speeds  (wind 
speeds  below  the  cut-in  speed),  analytical  models  which  accurately  charac- 
terize intermediate  and  high  speed  subintervals  are  most  desirable. 
Finally,  the  t  test  was  performed  again  to  investigate  whether 
there  is  a  significant  difference  between  the  mean  values  of  the 
power  ratios  of  the  beta  and  Weibull  distributions,  whose  parameters 
were  obtained  by  the  matching-moments  technique.   In  this  case  since 
the  means  of  two  distributions  are  compared  against  each  other,  rather 
than  comparing  one  mean  with  a  predetermined  value,  a  slightly  different 
form  of  Eq.  (2.6-2)  is  used,  namely  [19] 
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(2.6-3) 


where 


|N1  +  N2 


■j  +  N2 


;12  W1  -   1  +  s22  (N2  -  1) 


(K1   -  1)  +  (N2  -  1) 


K  and  N-  are  the  sample  sizes  of  each  distribution  to  be  compared, 
and  s  and  s   are  the  standard  deviations  of  each  distribution.   The 
quantity  s  is  called  the  pooled  estimate   of  the  population 
standard  deviation.   The  quantity  s   is  then  just  the  difference  in 
the  standard  deviations  of  the  two  means. 

The  results  are  shown  in  Table  2.6-6.   The  t  values  are  all  much 
greater  than  the  t  value  at  a  significance  level  of  0.001  and  an  infinite 
number  of  degrees  of  freedom.   Consequently,  it  can  be  stated  that  the 
two  distributions  do  yield  significantly  different  power  ratio  values 
and  hence  the  two  distributions  are  different.   From  Table  2.6-4,  it  can 
be  seen  the  beta  distribution  usually  overestimates  the  available  power 
since  the  power  ratios  corresponding  to  this  distribution  are  greater 
than  unity  whereas  the  Weibull  distribution  using  matching-moments  esti- 
mators underestimates  the  available  power  for  all  wind  turbine  sizes 
considered. 

In  summary,  analysis  of  544  observed  wind  speed  distributions  by 
the  routine  CURVEFIT  shows  that  either  the  Weibull  or  beta  distributions 
obtained  by  the  matching-moments  technique  gives  an  excellent  fit  to 
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2 
observed  windspeed  data.   Although  most  x  values  for  both  distributions 

are  too  large  to  accept  either  distribution  as  a  good  fit  at  even  the 

0.005  significance  level,  they  are  far  better  representations  of  the 

data  than  the  Keibull  distributions  obtained  by  using  a  least  squares 

2 
estimation  of  the  parameters.   These  unexpectedly  large  values  of  x 

for  the  Weibull  and  beta  distributions,  whose  parameters  are  estimated 
by  the  matching-moments  technique,  results  from  discrepancies  in  the  fit 
in  the  very  low  speed  subintervals.   The  remaining  wind  speed  sub intervals 
are  very  accurately  fit.   This  conclusion  is  supported  further  by  near 
unity  results  for  the  power  ratio  test,  which  places  large  emphasis  on 
intermediate  and  high  speeds.   Hence,  for  the  speed  subintervals  of 
interest,  the  Weibuli  or  beta  distribution  derived  from  the  matching- 
moments  method  are  generally  very  representative  of  the  wind  speed 
data. 
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3.   WIND  TURBINE  GENERATOR  OPTIMIZATION 

3.1   Introduction 

An  optimum  wind  turbine  generator  system  (WTGS)  depends  upon  many 
factors,  such  as  its  intended  application,  wind  characteristics,  economic 
considerations,  aesthetic  aspects,  system  reliability,  and  practical  con- 
straints placed  upon  its  design,  location,  building  materials,  etc. 
Consequently,  the  term  optimum   can  assume  many  meanings.   An  optimum 
system  may  be  one  which  is  completely  autonomous,  i.e.,  capable  of  gener- 
ating the  demanded  power  without  the  need  of  any  back-up  system.   Still 
another  interpretation  may  require  that  the  amount  of  time  the  wind 
machine  is  down  for  repairs  is  to  be  minimized.   For  this  case  the  optimum 
system  would  be  one  which  makes  use  of  the  most  reliable  components  in  the 
manufacture  of  the  wind  system. 

In  this  time  of  increasing  energy  awareness,  wind  power  is  being 
looked  upon  to  generate  electricity  on  both  a  central  station  and  a  local, 
decentralized  level.   As  noted  in  Chapter  1,  some  technical  and  economical 
problems  still  plague  the  large,  central  station  units.   Small  wind  turbines, 
on  the  other  hand,  have  produced  power  successfully  for  many  years.   In 
contrast  to  both  fossil  and  nuclear  fuels,  wind  power  is  both  a  relatively 
clean  and  renewable  source  of  energy.   However,  for  wind  generated  electricity 
to  make  a  significant  impact  as  a  decentralized  power  source,  it  must  be 
shown  to  be  capable  of  producing  energy  competitively  with  more  convent- 
ional sources  of  power.   For  decentralized  production,  economic  optimiz- 
ation becomes  a  key  factor.   Like  all  economic  decisions,  the  desirability 
of  one  alternative  over  another  is  based  upon  the  alternative  that  saves 
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the  user  the  most  money.   For  a  WTGS  to  be  economically  viable  it  must  be 
demonstrated  that  it  can  save  the  user  money,  otherwise  it  is  wiser 
economically  to  continue  to  purchase  all  of  the  demanded  power  from  the 
utility.   To  compute  the  economic  savings  afforded  by  use  of  a  certain 
size  WTGS  in  a  particular  application,  the  power  output  from  the  wind 
system  (given  wind  speed  data  at  the  location  where  the  WTGS  will  be 
used)  must  be  matched  with  the  load  demands  of  the  application  to  see  how 
much  of  the  load  the  wind  system  can  supply.   The  remainder  of  the  demand 
load  must  be  provided  by  a  back-up  system  or  purchased  from  an  energy 
utility.   In  addition,  the  capital  cost  of  the  WTGS  and  its  amortization 
becomes  an  economic  factor  to  be  considered  in  such  an  analysis. 

Most  previous  economic  studies  on  the  use  of  a  WTGS  to  generate 
power  have  simply  examined  the  total  amount  of  energy  generated  with  the 
assumption  that  the  generated  power  can  always  be  used.   However,  for 
decentralized  applications  the  demand  power  requirements  will  have  to  be 
matched  to  the  power  production  of  the  WTGS  both  of  which  will  vary 
throughout  the  day  and  from  season  to  season.   Very  little  work  has  been 
done  to  date  on  examining  the  economics  of  matching  a  dedicated  WTGS  to 
a  given  demand  load.   Developmental  Planning  and  Research  Associates  (DPRA) 
of  Manhattan,  Kansas,  has  recently  performed  pioneering  work  in  this  area 
[12].   Their  study  entailed  the  development  of  a  methodology  to  determine 
the  national  impact  that  economically  optimum  sized  WTGSs  could  have  in 
various  agricultural  applications.   However,  this  investigation  did  not 
examine  the  sensitivity  of  their  economic  optimization  procedure  to  changes 
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in  the  model  parameters,  e.g.,  mean  wind  speed,  variations  of  wind  speeds 
about  this  mean,  and  variations  in  load  demand.   Consequently,  in  the 
second  phase  of  this  work,  the  methodology  for  economic  optimization  of  a 
WTGS  given  wind  and  load  distribution  is  examined  and  the  results  of 
sensitivity  studies  on  the  optimization  model  are  reported.   In  addition, 
this  methodology  is  used  for  a  realistic  example  to  determine  the  WTGS 
size  needed  in  order  to  realize  maximum  savings  over  the  purchase  of  all 
of  the  required  power.   Also,  the  effect  of  detailed  wind  speed  and  load 
information  on  the  optimum  WTGS  size  is  illustrated  by  this  realistic 
example. 

Two  major  assumptions  have  been  made  about  the  application  of  a 
WTGS  to  a  dedicated  load  in  this  study.   First,  the  WTGS  was  assumed  to 
generate  AC  electric  power  compatable  with  that  supplied  by  the  utility 
grid.   Furthermore,  the  WTGS  was  connected  to  the  electric  grid  in  such 
a  manner  that  the  demanded  electric  power  which  cannot  be  supplied  by 
the  WTGS  was  purchased  from  the  electric  utility.   This  implies  the  use 
of  an  interfacing  system  between  the  WTGS  and  the  power  grid  so  that  a 
blending  of  wind  generated  and  utility  power  can  occur.   This  system 
will  keep  the  WTGS  in  synchronization  with  the  utility  grid  as  well  as 
monitor  whether  power  is  to  be  taken  from  the  grid  or  dumped  onto  the 
grid  depending  upon  the  demanded  load  and  the  output  of  the  WTGS.   Con- 
sequently, power  supplied  from  the  WTGS  will  be  indistinguishable  to 
the  user  from  power  supplied  by  the  utility  grid.   Second,  there  was 
no  energy  storage  capability  associated  with  the  WTGS.   Hence,  whenever 
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the  wind  was  of  sufficient  strength  to  generate  power  above  that 
demanded  by  the  load,  the  excess  was  either  wasted  or  sold  to  the 
utility  or  another  customer. 

3.2  Optimization  Methodology 

The  methodology  used  in  the  present  optimization  study  involves  the 
development  of  a  series  of  procedures  or  modules.   First,  the  power 
output  of  a  given  sized  WTGS  for  a  particular  distribution  of  wind 
speeds  must  be  matched  with  the  given  demand  load.   Because  demand 
loads  and  wind  speeds  typically  vary  throughout  a  day  and  throughout 
days  in  a  year,  load  demand  and  wind  speed  profiles  must  be  given  both 
for  various  times  throughout  a  day  and  for  many  typical  days  (or  seasons) 
throughout  a  year.   This  detail  in  the  wind  and  load  data  is  necessary  to 
compute  accurate  values  for  the  WTGS  performance  values,  e.g.,  the  generated, 
purchased,  and  excess  power.   Second,  once  these  WTGS  performance  values 
are  calculated,  the  savings  in  purchased  energy  costs  plus  the  profit 
from  selling  excess  power  (if  any)  must  be  compared  with  the  costs  assoc- 
iated with  the  installation  of  a  WTGS,  i.e.,  capital  and  maintenance  costs, 
to  find  if  the  given  WTGS  size  yields  a  net  savings.   Third,  other  WTGS 
sizes,  i.e.,  WTGSs  with  different  rated  powers  and  speeds,  must  be  evaluated 
in  order  to  find  the  WTGS  size  which  affords  the  maximum  economic  savings. 
If  too  small  a  WTGS  is  chosen  the  energy  production  will  be  insignificant 
compared  to  the  cost  of  the  WTGS.   Similarly,  if  too  large  a  system  is 
chosen,  the  WTGS  costs  will  overshadow  any  energy  savings.   Consequently, 
a  key  element  of  the  methodology  is  a  search  technique  to  optimize  the 
savings  achieved  by  installation  of  a  WTGS. 
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Because  of  the  inherent  structure  of  the  methodology  and  the  large 
number  of  different  sizes  of  VTGSs  which  will  have  to  be  tried,  this 
study  lends  itself  to  a  computer  algorithm.   A  flow  diagram  for  such  a 
computer  algorithm  is  shown  in  Fig.  3.2-1.   The  algorithm  consists  of 
combining  the  results  of  several  modules,  each  module  models  some  aspect 
of  necessary  information  needed  for  the  computation  of  the  size  parameters 
for  the  optimal  WTGS. 

For  the  computer  algorithm  to  give  accurate  estimates  of  the  WTGS 
needed  to  supply  a  particular  load  with  electrical  power,  the  input  data 
must  be  sufficiently  detailed.   This  includes  specifying  the  wind  speeds 
and  load  demand  requirements  for  several  time  intervals  throughout  each 
day  and  for  typical  days  throughout  a  year.   The  "typical"  days  used 
throughout  a  year  will  be  termed  seasons.      Within  a  given  season,  wind 
speeds  will  fluctuate  from  day  to  day  for  a  particular  daily  time  interval, 
and  representation  of  the  wind  speed  for  this  interval  by  a  constant 
value  for  every  day  in  the  season  is  not  at  all  realistic.   The  wind 
speed  in  any  daily  time  interval  must  be  characterized  by  a  wind  speed 
distribution.   Such  wind  speed  distribution  information  is  available  for 
many  locations  throughout  the  United  States  '[11,12].   Although  load  de- 
mands will  also  fluctuate  within  a  time  interval,  detailed  information 
about  these  fluctuations  is  generally  unavailable.   Consequently,  the  load 
demand  was  assumed  constant  for  every  time  interval.   This  assumption  of 
constant  load  demands  within  a  time  interval  is  not  altogether  inaccurate 
because  for  residences  with  an  established  living  pattern  or  an  enterprise 
with  well-defined  energy  applications,  the  load  will  vary  only  slightly 
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Fig.  3.2-1.   Block  diagram  of  optimization  methodology. 
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from  day  to  day  in  a  season  for  any  daily  time  period.   Finally,  the 
purchase  and  selling  price  (if  any)  of  electricity  must  be  established. 

With  the  optimisation  methodology  and  the  type  of  input  data 
specified,  models  for  each  of  the  program  modules  must  be  constructed. 
In  the  following  sections,  the  models  for  each  program  module  used  in 
this  investigation  are  described  in  detail. 

3.2.1  Wind  Models  and  Data  Requirements 

The  most  fundamental  wind  data  needed  for  this  optimization  study 
are  the  wind  speed  distributions  expected  for  a  given  daily  time  interval 
and  for  a  given  season.   If  detailed  wind  speed  data  are  available,  they 
can  be  used  directly  as  a  discrete  form  of  the  probability  density 
function  or  modeled  by  an  analytical  distribution.   Because  the  analytical 
wind  speed  distribution  is  conceptually  more  appealing,  since  in  reality 
the  wind  behaves  in  a  smoothly  varying  manner,  and  is  mathematically 
easier  to  handle,  the  Weibull  and  beta  distribution  functions  were 
examined  in  Chapter  2.   From  this  study  it  was  found  that  if  the  two 
parameters  of  either  of  these  analytical  representations  are  estimated 
by  a  matching-moments  technique,  accurate  characterizations  of  actual 
wind  speed  distributions  can  be  obtained.   The  use  of  analytical  distri- 
butions to  represent  the  wind  data  greatly  reduces  the  amount  of  data 
required  for  each  daily  time  interval  and  for  each  season  since  only  the 
two  parameters  of  the  distributions  need  be  specified. 

Since  many  wind  speed  distributions  are  needed  to  characterize  the 
variability  of  the  wind,  the  following  terminology  convention  is  adopted 
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for  the  remainder  of  this  report.  The  term  f   (v)  denotes  the  wind 
speed  distribution,  either  discrete  or  continuous,  in  the  i-th  daily 
time  interval  of  the  j-th  season.   The  number  of  time  intervals  and 
seasons  chosen  is  arbitrary  and  depends  upon  the  variability  of  the 
wind  throughout  the  day  and  the  year  at  the  site  chosen  for  the  WTGS, 
i.e.,  one  should  use  many  intervals  if  the  wind  distribution  changes 
rapidly  throughout  the  day  or  throughout  the  year. 

Often,  detailed  wind  speed  data,  usually  in  the  form  of  magnetic 
computer  tapes,  are  available  only  for  a  limited  number  of  locations. 
However,  the  required  wind  speed  distribution  functions  can  be 
synthesized  in  an  approximate  manner  from  less  detailed  meteorological 
data  than  are  readily  available.   For  such  a  synthesis,  two  pieces  of 
wind  data  are  required  for  each  location,  namely  (i)  the  wind  speed 
frequency  distribution  averaged  over  each  month  or  season,  and  (ii)  an 
average  wind  speed  distribution  as  a  function  of  time  of  day  for  each 
month  or  season.   The  overall  distribution  of  frequency  of  wind  speeds 
at  any  time  of  day  is  determined  primarily  by  the  local  weather  patterns 
of  fronts  and  other  slowly  varying  meteorological  phenomena.   Thus,  the 
relatively  short-time  diurnal  variation,  which  are  caused  primarily  by 
solar  heating  effects,  can  be  expected  to  affect  the  average  speed  at 
any  time  of  day  but  not  the  overall  shape  of  the  wind  frequency  distri- 
bution for  any  daily  time  interval.   Hence,  an  approximate  method  which 
can  be  used  to  generate  the  required  distribution  f  -1  (v)  is  to  assume 
its  overall  shape   (i.e.,  variance)  is  the  same  as  the  distribution  of 
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wind  speeds  in  the  j-th  season  but  with  the  mean  shifted  so  that  the 
mean  speed  corresponds  to  the  mean  speed  observed  in  the  i-th  time 
interval  of  the  j-th  season.   Consequently,  only  two  parameters  are 
needed  to  compute  the  necessary  wind  speed  distributions:   (i)  the  mean 
speed  in  the  i-th  time  interval  of  the  j-th  season,  and  (ii)  the  seasonal 
variance  of  wind  speeds  about  the  seasonal  mean  wind  speed.   Hence, 
with  the  mean  speed  and  variance  determined  for  each  daily  time  interval, 
the  matching-moments  technique  can  be  used  to  obtain  the  parameters  of 
an  analytical  distribution  for  each  daily  time  interval.   This  synthesis 
technique,  while  quite  straightforward,  is  only  approximate  and  should 
not  be  used  if  wind  speed  distribution  data  are  available  for  each  daily 
time  interval  for  each  season. 


3.2.2  WTGS  Response  Model 

The  response  function  for  the  WTGS  used  in  this  optimization  study 
is  an  idealized  one  which  was  discussed  in  Chapter  2,  namely 
0  ,  v  <  v  , 


R(v)  = 


P  (~)3 
r  v 


(2.5-8) 


"furl" 


As  noted  earlier,  v_   -  is  set  equal  to  v    ,  which  is  the  maximum 
furl  max 

speed  included  in  the  wind  speed  distribution  with  storm  conditions 
excluded.   Also,  cut-in  speed,  v  is  set  equal  to  0.46416  v  .   This 
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assumes  that  the  WTGS  is  idle  when  the  output  power  is  less  than  10% 
of  the  rated  power  for  the  WTGS. 

The  response  function,  R(v),  gives  the  net  electrical  power  output 
of  the  WTGS  exposed  to  a  steady  wind  speed  v.   Hence,  all  efficiency  and 
performance  factors  of  the  WTGS  are  taken  into  account  by  this  function  and 
no  explicit  description  is  needed  for  the  generic  type  of  wind  turbine 
or  electrical  conversion  system  used,  e.g.,  blade  diameter,  number  of 
blades,  aerodynamic  efficiency,  drive  train  efficiency,  and  inverter 
efficiency.   To  describe  completely  the  response  of  the  WTGS,  only  a 
pair  of  parameters,  P  and  v  ,  are  needed  to  describe  the  entire 
system  response. 

No  allowance  is  made  for  the  separate  effects  of  wind  gusts  on 
the  WTGS.   In  effect  it  is  assumed  that  the  wind  speed  distributions, 
f   (v),  have  been  determined  from  measurements  made  with  instruments 
which  have  the  same  dynamic  response  time  constants  as  the  idealized 
WTGS. 

In  most  realistic  systems,  there  is  also  a  furling  speed  above 

which  the  WTGS  is  feathered  or  shutdown  to  prevent  damage  from  high 

winds.   However,  for  most  wind  speed  distribution  data,  abnormally  high 

wind  conditions  are  suppressed.   By  suppression  of  such  storm  level 

winds  from  wind  speed  data,  the  existence  of  a  furling  speed  has 

effectively  been  incorporated  by  setting  f   (v)  to  zero  for  v  greater  than 

v    .   Hence,  v    becomes  the  furling  speed  of  the  model  WTGS  response 
max  max  or  r 

function. 
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The  assumed  response  function  used  in  this  study  resembles  the 
basic  features  of  many  response  functions  of  actual  wind  systems. 
However,  the  methodology  developed  in  this  study  is  sufficiently 
general  that  the  response  function  for  any  particular  WTGS  could  be 
readily  substituted. 

3.2.3  Matching  Load  to  Available  Power 

Once  the  wind  speed  model  and  WTGS  response  model  are  specified, 
the  output  power  that  a  particular  size  WTGS  can  generate  from  the 
given  wind  speed  conditions  can  be  matched  with  the  demand  load.   The 
three  power  calculations  needed  in  the  optimization  analysis,  are 
the  annual  required  electric  utility  energy  purchased,  the  annual  useful 
WTGS  energy  generated,  and  the  annual  excess  WTGS  energy  generated.   In 
the  following  subsections  the  formulas  used  to  calculate  each  of  the 
above  power  calculations  are  derived. 

(a)     Purchased  Electrical  Energy 

If  the  WTGS  is  too  small  or  the  wind  not  of  sufficient  strength 
to  handle  the  load  power  requirements,  electrical  power  must  be  purchased 
from  the  utility  to  supplement  the  power  generated  by  the  WTGS.   Hence, 
the  amount  of  power  that  must  be  purchased  over  the  i-th  daily  time 
interval  in  the  j-th  season  is 

P^  "  [P^3  -  R(v)]  H[P*3  -  R(v)],  (3.2-1) 

where 
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1,  x>0, 

0,  x<0. 


P  -1  is  the  electrical  power  demand  for  the  i-th  daily  time  interval 
d 

in  the  j-th  season,  and  R(v)  is  the  WTGS  electrical  power  output  when 
the  wind  speed  is  v.   The  unit  step  function,  H[x],  is  needed  to  describe 
the  case  when  the  power  demand  is  less  than  the  WTGS  power  generated, 
i.e.,  when  no  purchased  power  is  needed. 

Because  the  output  response  of  the  WTGS  depends  upon  the  wind 
frequency  distribution,  the  power  that  must  be  purchased  from  the  utility 
is  averaged  with  the  wind  frequency  distribution  to  compute  the  average 
or  expected  purchased  power  for  the  i-th  daily  time  interval  in  the  j-th 
season.   This  purchased  power  is  given  by 


vl3   =  P  [Py1  -  R(v)]  f1J(v)  H[P^  -  R(v)]dv. 
b    0   d  d 


(3.2-2) 


Since  high  wind  speeds  are  suppressed,  the  highest  wind  speed  observed 

at  a  location,  v    ,  can  replace  the  upper  limit  in  the  above  equation, 
max 

Also,  since  R(v)  is  equal  to  zero  when  v  is  less  than  v  ,  Eq.  (3.2-2) 
becomes 


vxJ   = 


I  max     ij    fiJ(v)dv  .  |  max  fiJ(v)dv,      P      <  M 

0  d  v  r~     ° 


/  P*j    flj(v)dv  -   /     R(v)flj(v)dv,      P      >  P^ 


(3.2-3a) 


(3.2-3b) 


where 


6   = 


mill    (v,,v        ),      v      <v,, 
d     max  c  —     d 
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z   =  min  (v   v   ) , 
d   max 


P   is  the  rated  power  of  the  WTGS,  and  v,  (which  is  less  than  v  )  is 
r  u  r 

the  speed  at  which  the  WTGS  first  produces  the  demanded  power.   The 
first  term  on  the  right  hand  side  in  Eq.  (3.2-3a)  can  be  simplified  by 
noting 


i   max  ij  .  . 
J      f   (v)dv 

0 


=  1  for  beta  distributions  and  actual 

wind  histograms, 
*  1  for  Weibull  distributions. 


(3.2-4) 


Therefore, 


r   max  _ij  ..ij  ,  .  , 
J     F  ,J  f  J  (v)dv 


ij  j  max  fij(v)dv  .  pij. 


(3.2-5) 


The  upper  limits  on  the  integrals  of  Eq.  (3.2-3b)  are  different 
than  those  of  Eq.  (3.2-3a)  because  when  P  is  greater  than  P   ,  there 
exists  a  demand  speed,  v  ,  such  that  P(v.)  is  equal  to  fA   The 
speed  v  is  calculated  by  solving  Eq.  (2.5-8)  for  v  ,  with  P(v.)  equal 


to  p: 


namely 


?1 

P 

r 


1/3 


(3.2-6) 


Two  sets  of  limits  are  needed  in  the  first  integral  of  Eq.  (3.2-3b) 

to  account  for  the  two  possible  relationships  between  v  and  v  .   When 

c      d 

v  is  greater  than  v.,  the  wind  speed  is  sufficiently  high  (according 
c  d 

to  Eq.  (3.2-6))  to  cause  the  WTGS  to  produce  power  but  is  still 
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below  the  cut-in  speed  and  the  WTGS  remains  idle.   Thus  power  must 

be  purchased  until  v  is  greater  than  v  .   When  v  is  less  than  or 
v  ■  c         c 

equal  to  v.,  power  is  only  purchased  until  the  wind  speed  reaches 

v.  and  the  WTGS  produces  enough  power  to  satisfy  the  demand  load, 
d 

For  wind  speeds  greater  than  v  ,  the  WTGS  generates  excess  power  and 

hence  no  power  needs  to  be  purchased.   If  in  any  of  the  above  integrals 

the  lower  limit  is  greater  than  the  upper  limit,  the  entire  Integral 
is  set  to  zero. 

Finally,  subsitutlon  of  the  explicit  form  for  P(v),  Eq.  (2.5-8), 
into  Eqs.  (3. 2- 3a)  and  (3.2-3b)  yields 


?' 


M 


-P  /  (^-)3  flj(v)dv-  P  f""fiJWdv,r  <  P^    (3.2-7.) 

r  ;   v  t  '  r  —  d 


/  P*j  f±j(v)dv  -  F  /   (^-)3  f   (v)dv,  Pr  >  P*j 

0  v    r 


(3.2-7b) 


Y  =  min  (v  ,  v   )  . 
r  max 


The  average  power  output  of  the  WTGS,  defined  by 


pU  =  / 


max  _,  .  ,ij  ,  ..  , 

R(v)  f   (v)dv, 


(3.2-8) 


cannot  always  be  used  to  estimate  the  expected  power  purchased  which 
intuitively  one  might  expect  to  be  simply  the  difference  between  the 


demand  and  average  WTGS  output  power,  i.e.,   P  -1  =  P  -1 


Pi;i.   Such  a 
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procedure  may  produce  an  erroneous  result  for  the  amount  of  purchased 
power.   Substitution  of  Eq.  (3.2-8)  for  the  second  integral  in  Eq.  (3.2-3a) 
yields 

pfj  =  P*j  -  Pij  ,  P  <  P^j  (3.2-9) 

b     d  r  —  d 

which  is  the  intuitive  result.   However,  if  Eq.  (3.2-8)  is  used  in 
Eq.  (3.2-3b),  the  result  is 


?t3   ■  P^J  -  P1J  +  f  [P(v)  -  P^J]  f1J(v)dv,  (3.2-10) 

r  a 


where 


t   =  max(v  .  ,v   ) . 

d  max 


Thus,  for  the  case  v  greater  than  v., 
*  r  to  d 

pij  >  pij  .  p-ij#  (3.2-11) 

b     d 

Hence,  the  purchased  power  would  be  underestimated  when  only  the 
average  power  output  of  the  WTGS  and  the  demanded  power  are  used. 

Finally,  because  utility  costs  are  based  on  energy  consumption,  the 
energy  used  in  every  time  interval  is  found  by  multiplying  the  expected 
power  purchased  in  that  interval  by  the  duration  of  the  time  interval, 
At . .   Then  this  quantity  is  summed  over  all  time  intervals  in  the  typical 
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day  of  a  season.   The  total  annual  purchased  energy  is  found  by  multiplying 
the  daily  consumption  by  the  number  of  days  in  a  season,  d.,  and  summing 
over  all  seasons.   Hence,  the  expected  total  annual  purchased  energy, 
E,  ,  is  given  by 

s      t  _.  . 
E =  I     d   I     ?}J    At  (3.2-12) 

b   j-1  J  i=l   b    X 

where  s  is  the  number  of  seasons  in  a  year  and  t  is  the  number  of  daily 
intervals. 

Cb)     Generated  Electrical  Energy 

The  annual  amount  of  electrical  energy  generated  by  the  WTGS  and  used 
by  the  load,E  ,  can  be  computed  from  the  result  of  the  previous  section  as 
the  difference  between  the  total  demand  and  total  purchased  electrical 
energy,  i.e., 


t 
E_  -  I     dj  I     P^3  At,  -  E,..         (3.2-13) 
j* 


=  I     dj  I     PdJ  Ati  -  V 
1=1  J  1=1 


(a)       Surplus  Electrical  Energy 

Whenever  the  WTGS  produces  more  power  than  is  demanded,  the  excess 
must  be  stored,  sold,  or  wasted.   Storage  has  not  been  considered  in 
this  study.   For  the  i-th  daily  time  interval  in  the  j-th  season  the 
amount  of  surplus  power  is 

Fgj  =  [R(v)  -  P*j]  H[R(v)  -  P^j].       (3.2-14) 
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Again,  the  unit  step  function  is  needed  so  when  R(v)  is  less  than  P,  , 
no  surplus  power  is  generated. 

Averaging  this  quantity  with  the  wind  speed  probability  density 
function  for  the  i-th  daily  time  interval  in  the  j-th  season  gives  the 
expected  surplus  power  (upon  substitution  for  R(v)) 


P   /  (— )3  flj(v)dv  +  P   /VmaX  fij(v)dv  -  pM  /Vmax  fij(v)  dv,  P  >  P*j 
r     v  r  a  r    a 

n    r  v  n 


?o  = 


C3.2-15a) 

,  P  <  P^  (3.2-15b) 

r  —  a 

where   n  =  max(v.,v  ), 
a   c 

and     G  =  min  (v  ,v   ) . 
r  max 

The  lower  limit,  max(v  ,v  ),  of  two  integrals  in  Eq.  (3.2-15a)  is 

required  since  there  can  be  no  surplus  power  generated  until  v  has 

been  reached. 

As  was  true  in  the  calculation  of  the  amount  of  purchased  electricity, 

the  use  of  only  an  average  power  output  from  the  WTGS  and  the  demand  power 

can  give  misleading  results  for  the  average  surplus  power.   This  derivation 

is  as  follows 

. .    v 

pi j  m   j  max  fiJ(v)[R(v)  _  p^]dv  (3.2-16a) 


=  /  maX  fij(v)[R(v)  -  P^jdv  -  /  f±j(v)  [R(v)  -  P^j]dv 
0 

(3.2-16b) 

_. .    . .   n  jj 

=  P1J  -  P*J  +  /   f1J(v)[P^J  -  R(v)]dv         (3.2-16c) 
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>  P1J  -  P*J  .  (3.2-16d) 

Hence,  the  amount  of  surplus  power  which  is  calculated  by  taking  the 
difference  between  the  average  WTGS  output  and  the  demand  power  will 
always  be  underestimated. 

The  annual  amount  of  excess  WTGS  energy,  E  ,  is  given  by 

s      t  _.  . 
E  =  "  d .   T  P13  At..  (3.2-17) 

j*l  J  i=l 

3.2.4  Evaluation  of  Integrals  for  P?"3  and  P13 
_ b s_ 

When  a  beta  or  Weibull  probability  density  function  is  substituted 
in  Eqs.  (3.2-7a),  (3.2-7b),  and  (3.2-15a)  the  integrals  cannot  be  evaluated 
analytically;  a  numerical  technique  must  be  used.   The  method  chosen 
in  this  study  is  Gauss-Legendre  quadrature.   The  quadrature  formula 
for  a  range  (a,b)  may  be  written  as  [20] 

b  r     _  -,   L      [-z.  (b  -  a)  +  b  +  a  \ 

j     g(x)dx  .  p-^-i   y  I  jM -2 j  ,  (3.2-18) 

a  ^     '    i=l     * 

where  L  is  the  number  of  quadrature  points  used,  w.  is  the  weighting 
factor  at  point  i,  and  z.    is  i-th  quadrature  ordinate.   Extensive 
tables  of  z.  and  w.  have  been  compiled  [21]. 

Because  the  beta  distribution  was  used  in  most  of  the  sensitivity 
studies  of  this  chapter  to  characterize  the  necessary  wind  speed 
distributions,  the  following  scheme  was  used  to  determine  a  quadrature 
order  that  would  minimize  computational  effort  yet  maintain  accuracy 
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of  integral  evaluation.   If  g(x)  is  a  polynomial  function,  then 
Gauss-Legendre  quadrature  will  evaluate  a  polynomial  of  the  order  2L  +  1 
or  less  exactly.   Since  the  beta  distribution  is  a  polynomial  function 

of  order  a  +  B  -  2,  if  o  and  B  are  integer  parameter  values  greater  than 

V2  i' 
unity,  integrals  of  the  form  /   f  J(v)dv  can  be  evaluated  exactly  by  a 

vl  vi  3  i' 

quadrature  order  of  !<(a  +  6  -   3)  and  integrals  of  the  form   /   (— -)   f  3  (v)dv 

Vl    r 
can  be  evaluated  by  a  quadrature  order  of  ^  (a  +  B) .   However, 

for  most  wind  speed  distributions,  the  parameters  of  the  beta  distribution, 
a  and  B,  are  non- integer  values,  but  by  choosing  the  quadrature  order 
to  be  the  nearest  integer  greater  than  h   (a  +  B) ,  it  can  be  expected  that 
the  numerical  evaluation  of  the  integrals  would  still  give  accurate, 
although  not  exact,  results.   However,  for  most  beta  distributions,  a 
quadrature  order  of  six  was  found  to  be  adequate  to  cover  most  values  of 
a  and  B  encountered  in  the  fitted  wind  distribution;  but,  before  using 
any  quadrature  order,  the  range  of  expected  a  and  B  values  should  be 
computed  and  if  large  values  are  indicated,  a  higher  quadrature  order 
should  be  used  for  those  cases  so  accurate  results  can  be  obtained. 

If  the  discrete  form  of  the  probability  density  function  is  used  in 
Eqs.  (3.2-7a),  (3.2-7b),  or  (3.2-15a),  the  integrals  can  be  evaluated 
exactly  since  the  discrete  distribution  is  a  constant  over  each  speed 
subinterval.   Consequently,  these  equations  must  be  divided  into  n 
speed  subintervals  and  integrated  between  the  endpoints  since  the  discrete 
distribution  assumes  a  different  but  constant  value  in  every  speed  sub- 
interval.   The  integrations  over  every  speed  subinterval  are  summed  to 
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vield  the  integration  over  the  entire  speed  range.   This  summation 
can  be  written  as  follows 

I      fjj   /  k+1  Q(v)dv,  (3.2-19) 

tc-1  k     v, 

k 

where  n  is  the  number  of  speed  intervals,  f  -1  (v)  is  the  value  at  the 
k-th  speed  subinterval  of  the  probability  density  function  for  the  i-th 
time  interval  in  the  j-th  season,  and  Q(v)  is  either  the  WTGS  response 
function,  R(v),  corresponding  to  the  speed  subinterval  to  be  integrated 
or  the  demand  power,  t      ,    for  the  i-th  time  interval  in  the  j-th  season. 
It  should  be  noted  that  Q(v)  is  a  polynomial  (or  constant)  that  can  be 
integrated  analytically. 

Therefore,  both  analytical  and  discrete  probability  density  functions 
can  be  used  to  evaluate  the  integrals  that  give  the  purchased  and  surplus 
electricity.   This  use  of  both  types  of  distributions  gives  flexibility 
to  the  optimization  program  in  that  the  effect  of  the  distribution  type 
on  the  optimum  WTGS  can  be  determined. 

3.2.5  WTGS  Cost  Model 

Because  WTGS  electrical  components  and  other  major  system 
components  have  not  been  produced  in  significant  quantitites,  the  costs 
of  various  sized  wind  turbines  are  difficult  to  determine.   The  costs 
vary  with  the  size  and  type  of  rotor,  gear  mechanism,  inverter,  and 
tower  as  well  as  the  type  of  feathering  device  used  to  prevent  damage 
to  the  WTGS  at  high  wind  speeds.   Furthermore,  control  units  which  auto- 
matically blend  WTGS  output  power  with  utility  power  to  meet  the  demand  load 
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have  not  been  commercially  developed.   A  compilation  of  costs  of  WTGSs 
already  in  production  has  been  made  by  DPRA  [12]  using  data  reported  by 
both  Obermeier  [10]  and  Rosen,  et.  al.  [22].   The  reported  values  are 
for  WTGSs  with  a  rated  speed  of  25  mph  and  are  shown  in  Fig.  3.2-2  in 
the  form  of  cost  per  kilowatt  ($/kW)  as  a  function  of  rated  power  (in 
kW) .   As  can  be  seen  from  Fig.  3.2-2,  the  data  from  Rosen,  et.  al. , 
yields  slightly  lower  capital  costs  than  the  data  from  Obermeier. 

To  interpolate  between  data  points  to  obtain  costs  for  various 
other  sizes  of  WTGSs,  a  quadratic  polynomial  has  been  fit  to  the  more 
optimistic  values  provided  by  Rosen,  et  al.   The  equation  for  this  cost 
model  is  [12] 

lnfel)  "  7.73971  -  0.46578[ln(P  )]  +  0.02573[ln(P  )2],      (3.2-20) 
KW  r  r 

F  >  1  kW.  v  ,   =  25  mph. 
r  —       ref 

where  v  .  is  the  reference  rated  speed   for  all  WTGSs  used  to  obtain 
ref 

this  result.   For  a  WTGS  rated  below  1  kW,  the  slope  of  Eq.  (3.2-20)  at 
P  equal  to  1  kW  is  used  to  give  the  simple  linear  cost  model 

4=  2297.8  (P  )-°-46578;  P  <  i  kW,  (3.2-21) 

kW  r  r 

v   ^  =  25  mph. 
ref 

Because  these  cost  models  assume  a  constant  rated  speed  (25  mph) 
and  because  the  optimization  methodology  varies  the  rated  speed  in 
order  to  find  the  optimum  combination  of  rated  power  and  rated  speed, 
a  conversion  factor  is  needed  so  that  the  values  obtained  from  the 
above  cost  models  can  be  applied  to  other  rated  speeds.   Eldridge  [1] 
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reports  that  capital  costs  of  a  WTGS  decrease  as  the  rated  speed  in- 
creases since,  for  a  given  output  capacity  of  a  WTGS,  the  size  and 
weight  of  the  blades,  rotor,  gear  train,  and  other  system  components 
will  generally  decrease  when  the  WTGS  is  designed  for  higher  rated 
speeds.   This  is  due  to  the  fact  that  power  derived  from  the  wind 
varies  linearly  with  the  area  swept  by  the  rotor  but  varies  with  the 
cube  of  the  wind  speed.   Hence,  for  a  given  rated  power,  area  decreases 
considerably  if  rated  speed  increases.   Using  data  from  references  10 
and  22,  DPRA  has  derived  an  expression  correcting  costs  (S/kW)  for 
different  rated  speeds.   For  a  WTGS  with  a  rated  speed  different  from 
that  used  to  derive  Eqs .  (3.2-20)  and  (3.2-21),  the  cost   can  be  corrected 
for  this  change  in  rated  speed  by  multiplying  the  capital  cost   ($/kW) 

by  (Wvre£)"2. 

Other  costs  associated  with  the  WTGS  are  those  of  operation  and 
maintenance.   Obermeier  [10]  reports  that  annual  operation  and  mainten- 
ance costs  are  a  fixed  percentage  of  the  capital  cost  of  the  WTGS  and 
depend  on  the  type  of  use  as  well  as  the  durability  of  system  components. 
Maintenance  costs  of  one  percent  for  factory  produced  wind  generators, 
three  percent  for  partly  owner  assembled  systems,  and  five  percent  for 
home  built  systems  are  typical.   In  this  study,  a  three  percent  annual 
maintenance  cost  was  assumed. 

Although  the  cost  model  described  above  and  used  in  this  study  is 
based  on  rather  uncertain  data  and  many  different  types  of  WTGS,  the 
model  is  felt  to  approximate  expected  costs  of  WTGS  given  a  sufficiently 
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large  degree  of  production.   However,  the  optimization  procedure  is 
written  such  that  any  other  cost  model  could  be  readily  substituted  for 
the  one  used  in  this  study. 

3.2.6   Calculation  of  Annual  Net  Savings 

The  net  financial  savings,  or  burden  if  money  stands  to  be  lost, 
realized  by  the  installation  of  a  WTGS  is  determined  by  adding  together  the 
money  saved  on  purchased  electrical  energy  costs  and  made  by  selling 
surplus  electricity  (if  any) ,  and  subtracting  from  this  value  the  costs 
associated  with  a  WTGS.   The  following  equation  expresses  this  relation- 
ship 

A  =  (E  )(C,  )  +  (E )(C J  -  AC,,T„  -  AOM,  (3.2-22) 

g    b       S    s       U'TGS 

where  A  is  the  annual  net  savings,  E  is  the  electrical  energy  (kWh) 

generated  in  a  year  by  the  WTGS  and  which  is  used,  E  is  the  surplus 

electrical  energy  generated  in  a  year  by  the  WTGS  which  is  wasted  or 

sold  to  another  user  (e.g.,  utility  or  other  enterprise),  C  is  the  cost 

($/kWh)  of  purchasing  electrical  energy  from  the  utility,  C  is  the  price 

for  which  surplus  WTGS  electricity  can  be  sold  (if  wasted,  C  equals 

zero),  ACTTrT^,„  is  the  effective  annual  cost  ($)  of  the  WTGS,  and  AOM 
WTGS 

is  the  effective  annual  cost  for  operating  the  WTGS. 

In  the  previous  section,  capital  cost  models  were  given  in  order 
to  compute  the  cost  of  installing  a  particular  size  WTGS.   Because  the 
WTGS  represents  such  a  large  investment,  the  capital  cost  is  amortized 
over  a  period  of  m  years,  at  a  yearly  interest  rate  i  .   The  fraction 
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of  the  initial  investment  required  each  year  for  the  amortization  is 
called  the  capital  recovery  factor   {CBS)    and  is  given  by  [23] 

i   (1  +  i  )m 

CRF  =  — —  .  (3.2-23) 

(1  +  ir)m  -  1 

The  capital  cost,  C  ,  of  the  WTGS  is  multiplied  by  the  CRF  to  obtain 
the  annual  net  cost,  A,  i.e. 

A  =  (C  )(CRF).  (3.2-24) 

c 

In  this  study,  the  interest  rate  was  assumed  to  be  10%  per  year 
with  a  20  year  WTGS  life  expectancy.   No  salvage  value  was  assumed  for 
the  WTGS  at  the  end  of  the  life  expectancy.   The  annual  operation  and 
maintenance  costs  for  the  WTGS  are  assumed  to  be  three  percent  of  the 
capital  cost  of  the  WTGS. 

3.2.7  Optimization  Technique 

Because  the  purpose  of  this  optimization  study  was  to  maximize 

the  annual  net  savings  described  by  Eq.  (3.2-22),  a  technique  to  find 

the  maximum  of  this  objective   function   was  needed.   The  net  annual 

savings  depends  upon  the  size  of  the  WTGS  through  its  rated  power  and 

rated  speed,  and  hence  the  problem  becomes  a  two-parameter  optimization 

problem.   An  algorithmic  technique  such  as  the  method  of  steepest  descent 

could  be  used  to  find  the  maximum  of  the  objective  function,  A(P  ,v  ). 

j  '    r'  r 

However ,  because  this  method  requires  evaluation  of  the  first  derivatives 
of  the  objective  function,  its  use  was  precluded  in  the  present  study 
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because  of  the  complexities  in  calculating  the  derivatives  of  the 
generated,  purchased,  and  surplus  powers.   Hence,  a  search  technique  is 
more  appropriate  since  such  a  method  is  merely  an  organized  procedure 
that  chooses  points  so  that  the  contours  of  the  two-dimensional  space 
are  scanned  in  order  to  find  the  point  that  yields  the  maximum  value 
of  the  objective  function.   The  technique  chosen  for  this  study  is  the 
sequential  simplex  pattern  search  method  [24,25,26].   The  computer  sub- 
routine used  to  calculate  the  maximum  by  the  simplex  technique  was  adapted 
from  a  program  written  by  Lai  [27]. 

3.3  Sensitivity  of  the  Optimally  Sized  WTGS  to  Problem  Parameters 

The  computer  code  BLOHARD  was  written  to  carry  out  the  optimization 
methodology  described  in  the  previous  sections  of  this  chapter.   This 
routine  is  based  on  a  similar  code  developed  by  DPRA  in  their  study  of 
wind  applications  in  agriculture  [12].   A  listing,  explanation,  and 
sample  output  of  routine  BLOHARD  are  contained  in  Appendix  C.   The  program 
BLOHARD  was  used  to  investigate  the  sensitivity  of  the  optimally  sized 
WTGS  to  various  problem  parameters  such  as  diurnal  variations  in  load, 
load  size,  diurnal  variations  in  wind  speed,  variance  in  wind  speed 
distributions,  variations  in  mean  wind  speed,  and  credit  given  for  surplus 
electricity.   Furthermore,  the  effects  of  seasonal  changes  in  the  wind 
speed  distributions  and  load  demands  are  also  examined. 

3.3.1  Model  Wind  and  Load  Profiles 

So  that  the  sensitivity  of  the  optimization  procedure  to  the 
problem  parameters  can  be  analyzed,  idealized  models  for  both  the  wind 
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speed  and  the  demand  load  must  be  developed.   Crawford,  et.al.  [28] 
report  that  due  to  solar  effects,  mean  wind  speeds  are  lowest  just 
after  midnight,  rise  gradually  after  sunrise  reaching  a  peak  in  the 
middle  afternoon,  and  then  fall  gradually  after  sunset  back  to  the 
low  value.   The  wind  speed  data,  which  were  analyzed  in  the  previous 
chapter,  and  were  obtained  from  the  National  Climatic  Center,  Asheville, 
North  Carolina,  for  17  locations  throughout  the  United  States,  exhibit 
such  diurnal  fluctuations.   One  simple  model  which  exhibits  such  diurnal 
variations  is  a  sinusoidal  variation  about  some  mean  daily  or  seasonal 
speed  of  the  form 

v(t)  ■  v  -  a  cos(-j7  -  r)>  (3.3-1) 

where  v(t)  is  the  average  or  expected  speed  at  time  t  (in  hours  beginning 
at  midnight),  v  is  the  mean  daily  speed,  and  "a"  is  the  amplitude  of  the 
variation  about  the  mean  speed.   Eq.  (3.3-1)  reaches  the  minimum  at  3  hours 
and  the  maximum  at  15  hours,  which  approximately  correspond  to  the  obser- 
vations made  by  Crawford,  et.al. 

To  generate  analytical  representations  of  such  model  wind  speed 
distributions,  e.g.,  wind  speed  distributions  for  eight,  three-hour 
intervals  in  the  day,  the  average  value  of  Eq.  (3.3-1)  was  calculated 
for  every  daily  interval  of  interest.   This  procedure  yields  an  average 
speed  for  every  daily  time  interval.   To  obtain  a  distribution  of  speeds 
about  this  mean  value,  for  each  daily  time  interval  a  dispersion  or 
variance  of  speeds  about  each  interval's  mean  speed  needs  to  be  specified. 
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For  the  model  winds  used  in  this  study  the  variance  was  assumed  to 
vary  linearly  with  the  interval's  mean  speed, 

s2  =  c  v.,  (3-3-2) 

l      l 

where  s.    is  the  wind  speed  variance  for  the  i-th  time  interval  of  the 
day,  V  is  the  average  speed  for  the  i-th  daily  time  interval,  and  c  is 
the  constant  of  proportionality  between  the  mean  speed  and  variance, 
called  the  coefficient  of  variation. 

Once  the  mean  and  variance  of  the  wind  speed  for  a  particular  daily 
time  period  are  calculated  from  Eqs.  (3.3-1)  and  (3.3-2),  a  model  wind 
speed  distribution  for  that  time  period  can  be  obtained  by  using  a  beta 
function  representation  whose  parameters  are  chosen  such  that  the  required 
mean  and  variance  are  realized.   In  Chapter  2  it  was  shown  that  such  a 
matching-moments  technique  was  both  analytically  very  simple  as  well  as 
vielding  beta  distributions  which  modeled  actual  wind  speed  data  very 
accurately.   To  obtain  model  wind  speed  profiles  for  use  in  the  sensitivity 
analyses,  each  day  was  divided  into  eight  time  intervals  or  periods  each  of 
three  hours  duration.   Diurnal  fluctuations  of  10%  and  25%  of  the  daily  mean 
wind  speed  were  chosen  (i.e.,  a  equal  to  0.1  v  and  0.25  v,  respectively  in 
Eq.  (3.3-1)).   Corotis  [7]  reports  that  diurnal  variations  of  about  10%  are 
representative  of  the  range  of  fluctuations  encountered  in  a  season  with  re- 
latively constant  winds  and  a  value  of  about  25%  is  typical  for  a  season 
with  rather  gusty  winds.   The  diurnal  fluctuations  as  a  function  of  daily 
time  interval  for  these  two  cases  are  listed  in  Table  3.3-1.   So  that  the 
effect  of  both  a  large  and  small  variance  in  the  average  speed  for  every 
daily  time  interval  can  be  studied,  coefficients  of  variation,  c,  equal  to 
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Table  3.3-1*  Average  Speed  as  a  Function  of  Time  Interval. 

Mean  Daily  Wind  Speed  is  10  Knots  for  Both  Cases. 


Time  Interval 
(hrs) 

Small  wind 

fluctuation 

10%  of  mean  speed 

(knots) 

Large  wind 
fluctuation 
25%  of  mean  speed 
(knots) 

0  - 

3 

9.10 

7.75 

3  - 

6 

9.10 

7.75 

6  - 

9 

9.63 

9.07 

9  - 

12 

10. A 

10.9 

12  - 

15 

10.9 

12.3 

15  - 

18 

10.9 

12.3 

18  - 

21 

10.4 

10.9 

21  - 

24 

9.63 

9.07 
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1  and  k    (corresponding  to  relatively  steady  and  gusty  wind  conditions, 

respectively)  were  selected.   Finally,  a  maximum  wind  speed  of  50  knots 

was  used  for  the  v    term  in  the  beta  distribution.   Normally,  speeds 
max 

above  50  knots  are  considered  storm  speeds  for  which  a  WTGS  would  not  be 
operated. 

Table  3.3-2  assigns  a  wind  model  number  to  the  different  permu- 
tations of  the  two  diurnal  variations  paired  with  the  two  coefficients 
of  variation  to  be  studied.   To  study  the  effect  of  seasonal  mean  wind 
speed  on  the  optimum  WTGS,  one  additional  case  was  run.   The  seasonal 
mean  speed  was  doubled  to  20  knots  and  combined  with  the  large  diurnal 
variation  and  large  coefficient  of  variation  to  form  a  model  characteristic 
of  high  speed,  gusty  wind  conditions. 

A  second  necessary  component  of  the  sensitivity  studies  is  the 
construction  of  model  demand  load  profiles.   DPRA  [12]  and  Obermeier  [10] 
report  that  typically  two  load  demand  peaks  are  reached  throughout  the 
day.   In  this  study,  the  greatest  load  demands  were  assumed  to  occur  between 
the  hours  of  9  and  12  and  again  from  hours  15  to  18.   Load  demands  were 
assumed  to  be  constant  throughout  every  three-hour  period,  i.e.,  there 
is  no  distribution  in  the  demand  about  the  given  load  for  a  particular 
time  interval.   Three  cases  were  considered  as  being  representative  of 
the  various  types  of  load  demands  that  may  occur  throughout  the  day; 
a  flashy   toad   (one  with  variations  in  the  mean  daily  load  of  100%) ,  a 
smooth   load   (one  with  variations  in  the  mean  daily  load  of  50%) ,  and  a 
constant   load   (no  variation  throughout  the  day).   Table  3.3-3  lists  the 
daily  variations  in  each  model  load  as  a  function  of  time  interval.   The 
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Table  3.3-2-  Wind  Model  Numbers  Assigned  to  the  Various 
Combinations  of  Diurnal  Fluctuations  and 
Coefficients  of  Variation. 


Wind  Model  Diurnal  Coefficient 

#  Fluctuation  of  Variation 


25%  of  seasonal 
mean  speed 

25%  of  seasonal 
mean  speed 

10%  of  seasonal 
mean  speed 

10%  of  seasonal 
mean  speed 

25%  of  seasonal 
mean  speed 
(20  knots) 
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Table  3.3-3,   Model  Load  Variations  as  a  Function  of  Time  Interval 
for  Flashy,  Smooth,  and  Constant  Loads.  (Ail  Loads 
are  normalized  to  a  daily  average  of  unit  demand  (1  kW)). 


Constant  Load 

1.0 
1.0 
1.0 
1.0 
1.0 
1.0 
1.0 
1.0 


Time  Interval 

Flashy  Load 

Smo 

oth  Load 

(hrs) 

1-3 

0.0 

0.5 

3  -  6 

0.5 

0.75 

6  -  9 

1.0 

1.0 

9-12 

2.0 

1.5 

12  -  15 

1.0 

1.0 

15  -  18 

2.0 

1.5 

18  -  21 

1.0 

1.0 

21  -  24 

0.5 

0.75 
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load  profiles  in  this  table  were  normalized  to  an  average  demand  of 
one  unit  (1  kW) ;  consequently,  for  a  given  mean  load,  the  values  in 
the  table  are  multiplied  by  the  desired  mean  daily  load  demand.   Mean 
daily  load  demands  of  5  kW  and  35  kW  were  chosen  as  values  representative 
of  a  small  residential  load  and  a  residential  load  combined  with  a  farm 
load  (or  possibly  a  combination  of  several  residential  loads),  respectively. 

3.3.2  Results  of  the  Model  Case  Studies 

The  model  wind  distributions  and  model  demand  load  profiles 
discussed  in  the  previous  section,  were  used  to  determine  the  sensitivity 
of  the  economically  optimum  WIGS  to  various  types  of  winds  and  loads.   In 
particular,  14  combinations  of  the  model  wind  speed  distributions  with 
model  loads  were  used.   These  combinations  or  cases   are  defined  in 
Table  3.3-4.   In  Cases  1  through  12  the  model  wind  profile  and  load  curve 
were  assumed  to  hold  for  the  entire  year,  i.e.,  only  a  single  season  was 
considered.   For  Case  13  two  seasons  were  used  for  the  year  while  Case  14 
used  four  seasons  with  actual  wind  speed  distributions  and  load  profiles 
characteristic  of  a  Kansas  farming  operation.   Although  the  single 
season  cases  are  unrealistic,  they  will  tend  to  accentuate  any  peculiar 
features  affecting  the  selection  of  the  optimum  WTGS.   Consequently, 
these  single  season  cases  are  useful  in  identifying  general  trends 
in  optimum  WTGS  size  without  specifying  large  amounts  of  wind  speed  or 
load  demand  data.   However,  the  interaction  of  seasonal  variations  in 
input  wind  speed  and  load  data  is  important  when  examining  the  compet- 
itiveness of  wind  energy  with  conventional  energy  sources;  hence,  more 
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Table  3.3-4. 

Case  Numbers  Assigned  to 

Various  Combinations  of  Wind 

and  Load  Models.   (S 

urpl* 

as  electrical  energy  is  assumed 

to  have  no  value  unl 

ess 

otherwise  noted) . 

Case  Number 

Wind  Model 

Type  of  Load 

(see  Table  3.3-2) 

(see  Table  3.3-3) 

1 

1 

Constant  -  5  kW  average 

2 

1 

Smooth  -  5  kW  average 

3 

1 

Flashy  -  5  kW  average 

4 

5 

Smooth  -  5  kW  average 

5 

1 

Smooth  -  5  kW  average 
Credit  for  surplus  (s  2c/kWh 

6 

4 

Flashy  -  35  kW  average 

7 

1 

Flashy  -  35  kW  average 

8 

3 

Flashy  -  35  kW  average 

9 

1 

Smooth  -  35  kW  average 

10 

1 

Smooth  -  5  kW  average 
Credit  for  surplus  @  2o/kWh 

11 

2 

Flashy  -  35  kW  average 

12 

5 

Smooth  -  35  kW  average 
Two  season  problem: 

13 

3  (season 

#1) 

Smooth  -  10  kW  average 

1  (season 

#2) 

Smooth  -  20  kW  average 

14 

Actual  seasonal   Four  Season  Problem: 

wind  data 

for 

Loads  typical  of  winter 

Wichita,  KS 

wheat/sorghum  farm 

operation 
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detailed  information  is  needed.   Cases  13  and  14  investigate  this  inter- 
action of  seasonal  variations  and  their  effect  upon  the  optimum  WTGS  size. 

A  summary  of  results  for  these  14  cases  is  presented  in  Tables 
3.3-5  through  3.3-7  for  various  assumed  costs  of  utility  supplied 
electricity.   The  results  in  these  tables  will  be  referred  to  throughout 
the  remainder  of  this  section  and  are  presented  here  for  convenience. 
The  optimum  WTGS  size,  i.e.,  rated  power  and  rated  speed,  is  listed  along 
with  the  amount  of  electricity  generated  and  used,  the  amount  of  electricity 
which  had  to  be  purchased  from  the  utility,  and  the  amount  of  surplus 
electricity  which  is  wasted  or  sold  to  other  users.   The  column  headed 
"%  self-sufficiency"  gives  the  percentage  of  the  entire  energy  demand 
that  is  generated  by  the  WTGS.   Finally,  the  annual  net  savings  of  the 
optimum  WTGS  is  listed.   In  the  following  subsections  the  results  of 
these  case  studies  are  discussed. 

(a)     Effect  of  Load  Variations  on  Optimum  WTGS 

Cases  1,2,3,7  and  9  best  illustrate  the  effect  of  the  three  different 
model  load  profiles  upon  the  optimum  WTGS.   Figs.  3.3-1  and  3.3-3  show 
that  the  optimum  rated  power  is  the  same  as  the  maximum  load  demand  if 
the  cost  of  utility  power  is  sufficiently  high.   However,  the  rated  speed 
increases  while  approaching  the  breakeven  value,  i.e.,  the  point  at 
which  costs  of  installing  a  WTGS  equal  the  fuel  costs  saved,  but  as 
purchased  electrical  costs  decrease  further  the  rated  power  may  begin  to 
decrease  rapidly.   This  behavior  is  expected  because  as  the  rated  speed  is 
increased,  initial  or  capital  costs  decrease  by  approximately  the  square  of 
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the  ratio  of  the  speed  change  due  to  the  decrease  in  the  VTGS  rotor 
size.   Capital  costs  are  also  reduced  if  the  rated  power  decreases, 
however,  in  the  5  kW  case,  the  saving  achieved  is  apparently  not 
significant  enough  to  compensate  for  the  loss  in  generation  capacity 
by  using  a  smaller  rated  power.   But  for  the  kK  mean  load  demand  cases, 
generator  rated  power  decreases  sharply  as  the  breakeven  value  is 
approached.   Initially,  as  the  cost  of  purchased  electricity  decreases, 
rated  speed  increases  to  decrease  capital  costs.   However,  an  upper  limit 
on  the  rated  speed  is  eventually  reached,  so  the  rated  power  decreases 
while  the  rated  speed  stays  relatively  constant.   As  electrical  costs  de- 
crease further,  the  WIGS  becomes  economically  infeasible  and  the  optimi- 
zation program  converges  to  a  "zero-cost"  size  WTGS,  i.e.,  one  with  a 
vanishingly  small  rated  power  with  an  exceedingly  large  rated  speed.   Such 
a  zero  cost  WTGS  is  one  that  is  physically  vanishingly  small. 

From  Figs.  3.3-2  and  3.3-4  it  can  be  seen  that  the  type  of  model 
load  has  a  relatively  small  effect  on  the  breakeven  cost  of  purchased 
electricity  (less  than  lc/kWh  difference).   At  all  costs  of  purchased 
electricity  that  yield  a  positive  annual  net  savings,  slightly  more 
savings  are  achieved  by  the  smooth  load  than  either  the  flashy  or  constant 
load.   Furthermore,  self-sufficiency  is  about  the  same  for  all  three  loads 
at  any  given  cost  of  purchased  electricity.   Consequently,  for  the  assumed 
diurnal  load  fluctuations,  only  a  small  increase  in  annual  net  savings 
and  decrease  in  the  breakeven  value  is  achieved  by  attempting  to  smooth 
the  variations  in  load  demand.   In  fact,   a  constant  load  fared  slightly 
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worse  than  a  smooth  load  as  far  as  annual  net  savings  and  self-sufficiency 
are  concerned.   This  intuitively  unexpected  result  is  a  consequence  of 
the  diurnal  variations  in  the  model  winds  used  for  these  examples  which 
are  inphase  with  the  variations  in  the  demand  load.   If  the  model  load 
fluctuations  occurred  in  the  time  intervals  with  the  lower  mean  speeds 
or  out-of-phase  with  the  model  wind,  the  constant  load  model  can  be 
expected  to  exhibit  a  larger  annual  net  savings  and  a  lower  breakeven 
value. 

One  peculiar  result  found  in  this  and  all  sensitivity  studies  in 
this  section  was  for  electrical  costs  slightly  below  breakeven,  an  optimum 
WTGS  was  found  by  BLOHARD  that  yielded  an  annual  net  savings  that  was  neg- 
ative, i.e,  a  loss.   This  result  is  caused  by  the  simplex  optimization 
technique  converging  upon  a  local  maximum  that  has  an  associated  negative 
savings  rather  than  converging  to  the  global  maximum  of  a  zero-size  WTGS 
with  an  annual  net  savings  of  zero.   However,  as  the  electrical  costs 
decrease  further  below  the  breakeven  cost,  the  optimization  routine 
converges  to  the  expected  annual  net  savings  of  zero.   Consequently,  such 
negative  results  for  net  savings  indicate  that  if  the  optimum  WTGS  must  be 
used  at  a  particular  electrical  cost  where  an  annual  net  loss  is  found,  this 
system  is  the  best  one  to  use,  since  this  system  will  lose  the  least 
amount  of  money. 

Furthermore,  the  simplistic  nature  of  the  curves,  i.e.,  horizontal 
or  vertical  trajectories,  shown  in  the  plots  of  rated  speed  and  rated 
power  as  a  function  of  cost  of  electricity,  is  a  consequence  of  using  a 
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single  season  to  describe  an  entire  year.   As  more  wind  speed  and  load 
data  are  supplied  so  that  seasonal  fluctuations  can  be  characterized, 
trajectories  are  obtained  that  exhibit  a  complex  variation  with  rated 
speed  and  power.   However,  sensitivity  analyses  using  a  single  season 
of  data  are  instructive  in  that  general  tendencies  in  the  optimum  WIGS 
size  can  be  identified. 

(b)     Effect  of  Average  Load  Size  on  the  Optimum  WTGS 

The  effect  of  average  load  size  on  the  optimum  WTGS  is  best  seen 
from  Figs.  3.3-2  and  3.3-4  which  show  annual  net  savings  as  a  function 
of  cost  of  purchased  electricity.   Larger  average  loads  have  approximately 
a  1.5  times  reduction  in  the  breakeven  electrical  cost  for  the  same  wind 
and  normalized  load  models.   Breakeven  costs  for  the  35  kW  average  load 
range  between  5.5c/kWh  to  6.5c/kWh  while  the  5  kW  average  loads  have 
breakeven  costs  of  about  10.5c/kWh.   This  variation  with  average  load 
is  expected  since  larger  WTGSs  are  less  expensive  on  a  per  unit  capacity 
basis  than  are  smaller  systems.   Consequently,  because  of  the  large  net 
savings  potential,  large  loads  are  more  attractive  to  wind  energy  appli- 
cations.  However,  as  the  WTGS  rated  power  increases,  a  value  of  rated 
power  is  reached  whereby  the  cost  function  assumed  in  Eq.  (3.2-20)  is 
no  longer  valid.   This  is  due  to  the  requirement  of  more  expensive  and 
sophisticated  control  systems  and  other  components  for  a  large  WTGS  that 
are  not  needed  for  a  small  WTGS.   Hence,  for  very  large  loads,  the  WTGS  cost 
function  should  be  modified  for  large  units  to  reflect  the  expected  in- 
crease in  per  unit  capacity  costs  with  size. 
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Besides  lower  breakeven  electrical  costs,  the  cases  with  large 
average  load  demands  exhibit  a  greater  self-sufficiency  at  a  specified 
electrical  cost  than  the  lower  average  load  demand,  given  the  same  wind 
models.   For  example,  Cases  2  and  9,  3  and  7,  and  5  and  10  are  the 
same  wind  models  and  normalized  load  models  with  the  only  difference 
lying  in  the  average  load  demand.   From  Table  3.3-6,  which  shows  the  re- 
sults of  the  optimum  WTGS  for  the  same  cost  of  purchased  electricity,  the 
self-sufficiency  factor  is  at  least  17%  greater  for  the  large  average  load 
demand  than  for  the  small  one.   Hence,  this  too  shows  larger  average 
loads  are  better  candidates  for  wind  energy  applications  than  small 
average  loads. 

As  the  breakeven  electrical  cost  is  approached,  the  rated  power 
for  the  cases  with  the  large  average  load  demand  is  very  sensitive  to 
slight  changes  in  electrical  costs.   From  this  behavior,  it  can  be  in- 
ferred that  large  loads  do  not  lose  a  significant  amount  of  generation 
capacity  by  reducing  the  generator  rated  power  in  order  to  reduce  capital 
costs. 

(a)     Effect  of  Diurnal  Variations  on  the  Optimum  WTGS 

Cases  7  and  8  demonstrate  the  effect  of  diurnal  variations  on 
the  optimum  WTGS.   From  the  results  shown  in  Fig.  3.3-6,  there  is  only 
a  slight  effect  on  annual  net  savings;  the  case  with  larger  diurnal 
variations  saves  more  money  than  the  case  with  small  diurnal  variations. 
However,  the  differences  in  annual  net  savings  become  negligible  with 
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decreasing  electrical  costs  and  both  cases  have  comparable  breakeven 
costs,  i.e.,  less  than  0.5c/kWh  difference.   Figure  3.3-5  indicates 
that  the  optimum  WTGS  size  is  the  same  for  both  cases  up  to  a  certain 
electrical  energy  cost.   Below  this  critical  cost,  the  case  with  the 
small  diurnal  variations  shows  sharp  reductions  in  the  rated  power  as 
electrical  costs  decrease  whereas  the  case  with  large  diurnal  variations 
remains  at  the  same  rated  power  and  just  increases  the  rated  speed.   As 
breakeven  is  approached,  both  systems  have  reduced  power  ratings,  but  the 
case  with  large  diurnal  variations  has  higher  rated  speeds. 

An  explanation  of  this  observation  is  that  there  is  less  variability 
in  the  wind  speeds  that  occur  during  the  daily  wind  speed  peak  in  the 
case  with  a  small  diurnal  variation  than  the  case  with  a  large  diurnal 
variation.   This  occurs  because  the  variance  of  wind  models  used  is 
assumed  to  be  directly  proportional  to  the  average  wind  speed  for  any 
daily  interval  (see  Eq .  (3.3-2)).   Since  the  daily  wind  speed  peak 
coincides  with  the  larger  load  demands  for  these  model  cases,  there  is 
a  greater  probability  of  higher  wind  speeds  occurring  in  the  daily 
wind  speed  peak  of  the  wind  model  with  the  larger  diurnal  variations. 
Hence,  because  of  this  greater  probability  of  higher  speeds,  the  rated 
speed  of  the  optimum  WTGS  can  be  increased  (so  as  to  lower  capital  costs), 
and  still  provide  a  sufficient  amount  of  electricity  to  try  to  cover  the 
peak  load  demand,  rather  than  lowering  capital  costs  by  decreasing  the 
rated  power  and  thereby  decreasing  the  capability  of  the  WTGS  to  supply 
the  peak  load  demand.   However,  this  probability  of  higher  than  average 
wind  speeds  does  not  occur  when  there  is  only  a  small  diurnal  variation; 
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hence,  the  rated  power  must  be  reduced  instead,  so  that  capital  costs  can 
be  lowered  and  an  economically  feasible  system  obtained. 

In  addition  to  having  only  a  slight  effect  on  the  annual  net 
savings,  diurnal  variations  affect  the  self-sufficiency  of  the  WTGS 
only  slightly.   Winds  with  large  diurnal  variations  have  a  somewhat 
larger  self-sufficiency  than  winds  with  small  diurnal  variations  about 
the  same  seasonal  mean  wind  speed.   This  effect  is  seen  by  examining 
the  results  of  Cases  7  and  8  listed  in  Tables  3.3-6  and  3.3-7.   Con- 
sequently, little  effect  on  the  optimum  WTGS  is  observed  by  diurnal 
variations   if  the  overall  mean  wind  speed  remains  constant. 

(d)     Effect  of  Wind  Speed  Fluctuations  on  the  Optimum  WTGS 

The  effects  of  changes  in  the  variance  of  the  wind  speed  distribution 
of  speeds  about  the  average  speed  in  a  particular  time  interval  has  also 
been  investigated.   In  the  model  wind  profiles,  coefficients  of 
variation  with  values  of  four  and  unity  are  used.   Cases  6  and  8  and 
Cases  7  and  11  are  excellent  examples  because  the  only  difference 
between  these  pairs  is  in  the  values  used  for  the  coefficients  of 
variation  in  the  wind  profiles.   Figure  3.3-7  shows  that  for  large 
electrical  costs,  the  annual  net  savings  for  the  cases  with  small  variances 
and  either  large  or  small  diurnal  variations  are  quite  large.   However,  as 
the  cost  of  electricity  decreases,  the  curves  drop  off  sharply,  ending 
up  with  breakeven  values  that  are  nearly  the  highest  of  all  cases  studied. 
Figure  3.3-8  shows  the  effects  differences  in  the  variances  of  wind 
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speed  distributions  have  on  the  size  of  the  optimum  WTGS .   As  expected, 
for  sufficiently  high  electrical  costs,  the  rated  power  attains  the  value 
of  the  largest  load  demanded  in  the  season.   The  cases  with  a  small  co- 
efficient of  variation  exhibit  much  lower  rated  speeds  (and  have  higher 
capital  costs).   Again,  with  decreasing  electrical  costs,  the  rated  speed 
increases  while  rated  power  remains  constant,  along  with  the  characteristic 
sharp  break  as  rated  speed  stays  constant  while  rated  power  decreases 
rapidly.   However,  the  range  of  rated  speeds  at  which  the  rated  power 
remains  constant  is  significantly  smaller  in  the  cases  with  the  small 
coefficient  of  variation  than  the  cases  with  the  large  coefficient  of 
variation. 

These  differences  between  the  cases  caused  by  different  wind  speed 
variances  can  best  be  explained  by  noting  that  the  wind  speed  distribution 
with  a  coefficient  of  variation  of  unity  is  very  highly  peaked  around  the 
average  value.   Because  of  this  small  variation  in  wind  speeds  about  the 
mean,  there  is  a  relatively  small  probability  of  wind  speeds  much  greater 
than  one  standard  deviation  beyond  the  mean  speed.   Hence,  although  the 
optimum  rated  speed  increases  with  decreasing  cost  of  purchased  electricity, 
too  much  of  an  increase  in  rated  speed  would  place  the  rated  speed  value 
at  a  point  on  the  tail  of  the  wind  speed  distribution  where  the  wind  has  only 
a  slight  probability  of  occurring.   Since  very  little  power  could  be 
generated  by  increasing  the  rated  speed,  the  rated  power  must  be  reduced 
instead  in  order  to  lower  the  capital  costs  of  the  WTGS.  Therefore, 
the  WTGS  becomes  infeasible  at  a  much  lower  rated  speed  for  a 
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wind  speed  distribution  with  a  small  coefficient  of  variation  than  for 
a  distribution  with  a  large  coefficient  of  variation  because  the  system's 
rated  capacity  must  be  reduced  sooner  when  the  wind  speed  distribution 
has  a  smaller  coefficient  of  variation. 

Although  the  breakeven  costs  are  higher  for  a  wind  model  with  a 
small  coefficient  of  variation,  the  self-sufficiency  is  enhanced.   Tables 
3.3-6  and  3.3-7  show  that  the  cases  under  study  have  values  of  self- 
sufficiency  that  differ  from  each  other  by  at  least  15%  for  equivalent 
costs  of  electricity.   Hence,  although  breakeven  costs  are  higher,  wind 
speeds  with  distributions  that  are  grouped  very  closely  about  the  mean 
speed  can  produce  much  greater  amounts  of  power  than  wind  speed  distributions 
that  are  highly  dispersed  about  the  same  mean  speed. 

(e)     Effect  of  Doubling  Mean  Wind  Speed  on  the  Optimum  WTGS 

To  study  the  effect  of  variations  in  the  seasonal  mean  speed  upon 
the  optimum  WTGS,  a  model  wind  speed  was  doubled  from  10  knots  to  20  knots. 
The  wind  model  used  had  both  the  large  diurnal  variations  and  the  large 
coefficient  of  variation,  while  the  model  load  chosen  was  the  flashy  load 
with  average  load  demands  of  both  5  kW  and  35  kW. 

The  results  for  Cases  4  and  12  (shown  in  Fig.  3.3-9)  indicate 
that  annual  net  savings  and  breakeven  costs  are  dramatically  increased 
and  decreased,  respectively.   In  increasing  the  mean  wind  speed  from  10 
to  20  knots,  the  breakeven  value  for  the  5  kW  average  load  model  was 
reduced  to  about  2.75c/kWh  from  the  previous  value  of  10.5c/kWh  and  for 
the  35  kW  average  load  model,  the  breakeven  value  was  reduced  to  1.5?/kWh 
from  the  previous  value  of  6.0c/kWh.   Tables  3.3-5  through  3.3-7 
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show  that  self-sufficiency  was  increased  also.   Even  at  the  breakeven 
value,  the  optimum  WTGS  for  both  the  large  and  small  average  load  demand 
has  a  self-sufficiency  in  excess  of  35%.   Figure  3.3-10  shows  that  results 
similar  to  the  10  knot  seasonal  mean  speed  cases  (Cases  2  and  10)  are  ob- 
tained, i.e.,  rated  power  achieves  the  value  of  the  maximum  load  demand 
and  with  decreasing  electrical  costs  the  rated  speed  increases  until 
breakeven  is  reached  (like  the  previous  5  kW  mean  load  demand  cases) 
or  the  rated  speed  increases  to  a  certain  point  and  then  stays  constant 
while  the  rated  power  decreases  until  breakeven  is  reached  (similar  to  the 
previous  35  kW  mean  load  demand  cases).   Furthermore,  the  range  of  speeds 
at  which  the  rated  power  of  the  WTGS  remains  constant  is  much  greater  than 
the  corresponding  case  with  a  10  knot  mean  wind  speed.   This  is  to  be  ex- 
pected since  for  the  same  coefficient  of  variation,  a  mean  speed  of  20  knots 
produces  a  greater  variability  of  wind  speeds  about  this  mean.   Hence, 
because  there  is  a  greater  probability  of  higher  wind  speeds  occurring 
for  a  mean  wind  speed  of  20  knots  than  for  one  of  10  knots,  the  rated 
speed  of  the  optimum  WTGS  can  be  increased  to  lower  capital  costs  of  the 
WTGS  before  reaching  the  point  where  rated  power  must  be  decreased  to 
lower  capital  costs. 

The  investigation  of  the  effect  of  increases  in  seasonal  mean  speeds 
on  the  optimum  WTGS  is  important  because  seasonal  mean  speeds  generally 
increase  with  the  height  of  a  WTGS  above  ground  [3].   Much  more  power 
can  be  generated  by  increasing  the  height  of  the  WTGS  above  ground. 
However,  because  increasing  the  tower  height  increases  the  cost  of  the 
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Fig.  3.3-10.   Size  of  optimum  WTGS  as  a  function  of  cost  of  electricity 
Cases  2,  4,  9,  and  12. 
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WTGS,  there  will  be  a  tradeoff  in  increased  power  production  and 
increased  tower  capital  costs. 

(f)     Effect  of  Credit  for  Surplus  Electricity  on  the  Optimum  WTGS 

Of  considerable  interest  in  the  development  of  the  WTGS  is  the 
effect  of  the  utility  giving  credit  for  any  surplus  energy  that  is  fed 
back  into  the  utility  grid.   A  credit  of  2c/kWh  was  chosen  since  any 
credit  given  by  the  utility  in  the  near  future  can  be  expected  to  be 
considerably  less  than  the  price  that  is  charged  to  purchase  electricity. 
This  expected  small  credit  is  due  to  the  fact  that  costs  of  purchased 
electricity  reflect  the  capital  costs  of  the  utility,  amortization  of 
the  transmission  and  distribution  systems,  costs  of  protective  systems, 
etc.,    in  addition  to  actual  fuel  costs.   Since  any  power  fed  into  the 
utility  grid  by  the  WTGS  will  save  the  utility  only  the  cost  of  fuel, 
the  low  credit  value  is  reflective  of  this  partial  savings.   This  credit 
was  applied  to  both  the  5  kW  and  35  kW  average  demand  load  models. 

The  effects  of  credit  for  surplus  electricity  are  best  seen  by 
contrasting  Cases  2  and  5  and  Cases  9  and  10.   The  only  difference 
between  each  of  these  cases  is  that  Cases  5  and  10  receive  credit  for 
surplus  electricity.   The  results  obtained  for  the  optimum  WTGS  are 
tabulated  in  Tables  3.3-5  through  3.3-7.   As  can  be  seen,  considerably 
larger  optimal  rated  powers  are  attained  in  order  to  generate  more  sur- 
plus energy.   However,  the  rated  speeds  remain  about  the  same  as  in  the 
cases  where  no  credit  was  given.   Furthermore,  only  a  moderate  increase  in 
self-sufficiency  is  achieved  if  credit  is  given  for  surplus  electricity. 
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Figures  3.3-12  and  3.3-13  show  that  there  is  only  a  slight  increase 
in  the  annual  net  savings  but  very  little  effect  in  the  breakeven  cost. 
Finally,  as  electrical  costs  decrease,  Figs.  3.3-11  and  3.3-14  show 
that  both  rated  power  and  rated  speed  decrease  in  a  smooth  but  complex 
manner,  unlike  the  sharp  breaks  in  the  trajectories  of  rated  speed  versus 
rated  power  in  the  cases  where  no  credit  is  given  for  surplus  energy. 
Consequently,  for  the  expected  low  value  of  credit  received  for 
surplus  power,  the  effect  on  the  optimum  WTGS  is  slight.   Although  more 
net  savings  are  achieved  at  every  cost  of  purchased  electricity,  there 
are  only  slight  improvements  in  the  breakeven  values  and  self-sufficiencies. 
Unless  larger  credits  are  given  for  generation  of  surplus  energy  by  a 
WTGS,  it  is  conjectured  that  more  effective  use  of  this  surplus  energy 
can  be  made  by  storing  the  excess,  e.g.,  in  batteries. 

(g)     Effect  of  Seasonal  Variations  on  the  Optimum  WTGS 
(i)   Two  Seasons: 

In  the  first  multiple  season  case  study,  two  seasons  of  load  and 
wind  data  were  used  to  characterize  the  entire  year.   For  the  first 
season,  wind  model  3  supplies  power  to  a  smooth  load  with  a  10  kW 
seasonal  average  and  for  the  second  season,  wind  model  1  supplies 
power  to  a  smooth  load  with  a  20  kW  seasonal  average.   Figure  3.3-15 
shows  the  results  of  this  two  season  optimization  problem  along  with 
the  results  of  each  season  run  as  a  single  season  case,  i.e.,  used  for 
the  entire  year.   As  can  be  seen,  for  high  costs  of  purchased  electricity, 
the  rated  power  of  the  optimum  WTGS  is  the  same  as  the  maximum  load  demand 
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for  the  entire  year.   As  the  cost  of  electricity  decreases  the  optimal 
rated  speed  increases  while  holding  the  rated  power  constant  in  order  to 
decrease  capital  costs.   Again,  the  tendency  is  to  try  to  cover  as  much 
of  the  peak  demand  as  possible.   However,  if  this  cannot  be  achieved,  the 
rated  speed  remains  constant  while  rated  power  decreases  so  as  to 
decrease  capital  costs.   But  because  there  are  two  seasons  with  different 
average  load  demands,  there  exists  another  rated  power  at  which  the 
WTGS  can  sufficiently  supply  the  peak  of  the  lesser  load  demand  as 
well  as  supply  an  adequate  portion  of  the  greater  load  demand,  though 
it  cannot  cover  all  of  the  maximum  demand  load.   Hence,  as  the  cost  of 
electricity  approaches  breakeven,  the  size  of  the  WTGS  remains  at  this 
lower  value  of  rated  power  while  rated  speed  continues  to  increase. 

Figure  3.3-16  shows  that  the  annual  net  savings  for  the  two 
season  case  lies  between  the  values  obtained  when  each  season  characterizes 
an  entire  year.   Furthermore,  the  breakeven  cost  is  approximately  the 
average  value  of  the  individual  seasons  when  each  season  is  used  for 
the  entire  year.   Consequently,  the  effect  of  seasonal  variations  on 
the  optimum  net  savings  is  that  of  averaging  the  results  obtained  from 
the  optimization  procedure  when  each  season  is  run  separately. 

Although  Fig.  3.3-15  still  exhibits  the  step-like  characteristic 
of   the  single  season  models,  it  must  be  kept  in  mind  that  the  seasonal 
load  variation  also  exhibited  these  characteristics,  i.e.,  the  two 
seasonal  load  profiles  differ  by  a  factor  of  two.   In  realistic  cases, 
the  seasonal  wind  and  load  changes  are  not  as  pronounced,  but  vary 
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smoothly  from  season  to  season.   Hence,  as  more  seasons  are  added 
it  can  be  expected  that  the  trajectory  of  the  WTGS  size  as  a  function 
of  cost  of  electricity  will  vary  in  a  smooth  but  more  complex  manner, 
(ii)  Four  Seasons: 

A  four  season  case  was  studied  to  see  if  the  results  obtained  from 
the  idealized  model  wind  speed  and  load  models  could  be  used  to  predict 
the  behavior  of  realistic  wind  speed  and  load  profiles.   A  load  model 
representative  of  a  Kansas  winter  wheat-sorghum  farm  was  chosen.   Wind 
speed  data  for  Wichita,  Kansas  was  obtained  from  the  National  Climatic 
Center  and  are  listed  in  Table  3.3-8.   These  wind  speed  data  are  assumed 
to  be  applicable  to  winter  wheat-sorghum  operations  near  this  city.   Load 
information  for  the  residence  and  farm  was  obtained  from  DPRA  [12].   Both 
the  residence  and  farm  load  are  shown  in  Tables  3.3-9  and  3.3-10.   The 
residential  load  is  based  on  a  totally  electric  home  for  a  family  of  four 
under  a  normal  daily  living  pattern.   The  only  significant  farm  load 
comes  from  aeration  fans  which  are  used  in  the  wheat  storage  bins  to 
prevent  spoilage. 

The  discrete  wind  speed  distributions  in  every  three-hour  interval 
are  modeled  by  the  beta  distribution  using  the  method  of  matching- 
moments .   The  parameters  of  the  approximating  beta  distribution  of  the 
wind  distributions  for  each  daily  and  seasonal  time  interval  are  tabulated 
in  Table  3.3-11.   As  a  comparison  of  the  effect  of  using  approximating 
analytical  distributions,  the  discrete  wind  speed  distributions  were  also 
used  in  the  optimization  routine  for  this  realistic  case.   The  results  of 
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Table  3.3-9.   Typical  Residence  Load  for  Family  of  Four  on 
Western  Kansas  Winter  Wheat-Sorghum  Farm 
(From  Ref .  12). 


Time 

Interval 

(hrs.) 

Season 

0-3 

3-6 

6-9 

9-12 

12  -  15 

15  -  18 

18  -  21 

21  -  24 

Load 

Demand  (kW) 

Fall 

5.96 

3.77 

15.6 

8.48 

6.47 

7.53 

15.9 

12.6 

Winter 

19.7 

20.0 

31.8 

22.2 

20.2 

21.0 

36.9 

28.6 

Spring 

11.5 

12.0 

21.1 

11.2 

9.22 

4.78 

15.9 

12.6 

Fall 

7.91 

5.72 

19.9 

12.8 

17.8 

16.13 

24.5 

18.9 

Table  3.3-10.   Typical  Load  for  Kansas  Winter  Wheat-Sorghum  Farm 
(From  Ref.  12). 


Season 

0-3 

3-6 

6-9 

Time 
9-12 

Interval 

12  -  15 

(hrs.) 
15  -  18 

18  -  21 

21  -  24 

Load  Demand 

(kW) 

Fall 

4.5 

4.5 

4.5 

4.5 

4.5 

4.5 

4.5 

4.5 

Winter 

4.5 

4.5 

4.5 

4.5 

4.5 

4.5 

4.5 

4.5 

Spring 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

Fall 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 
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Table  3 

.3-11.   Parameters 
Wind  Speed 

of  Beta 

Data. 

Distribution  for 

Wichita,  Kansas 

Time 

V 
max 

Seaon 

Interval 

a 

B 

(hrs) 

(knots) 

Fall 

0 

-   3 

1.32 

3.73 

33.5 

3 

-  6 

1.56 

4.42 

33.5 

6 

-  9 

1.67 

4.38 

33.5 

9 

-  12 

2.14 

4.13 

33.5 

12 

-  15 

2.28 

5.21 

40.5 

15 

-  18 

2.00 

3.72 

33.5 

18 

-  21 

1.41 

2.94 

27.5 

21 

-  24 

1.53 

4.14 

33.5 

Winter 

0 

-  3 

1.59 

3.72 

33.5 

3 

-  6 

1.54 

3.75 

33.5 

6 

-  9 

1.56 

3.71 

33.5 

9 

-  12 

1.85 

3.64 

33.5 

12 

-  15 

2.00 

3.50 

33.5 

15 

-  18 

2.09 

3.49 

33.5 

18 

-  21 

1.57 

3.88 

33.5 

21 

-  24 

1.60 

3.74 

33.5 

Spring 

0 

-   3 

1.61 

3.34 

33.5 

3 

-  6 

1.99 

5.44 

40.5 

6 

-  9 

2.05 

3.88 

33.5 

9 

-  12 

2.68 

5.07 

40.5 

12 

-  15 

2.81 

6.34 

47.5 

15 

-  18 

2.45 

4.45 

40.5 

18 

-  21 

2.07 

5.20 

40.5 

21 

-  24 

1.55 

3.32 

33.5 

Summer 

0 

-  3 

1.88 

4.14 

27.5 

3 

-  6 

1.74 

4.46 

27.5 

6 

-  9 

1.90 

4.12 

27.5 

9 

-  12 

2.41 

5.49 

33.5 

12 

-  15 

2.36 

5.08 

33.5 

15 

-  18 

2.51 

5.20 

33.5 

18 

-  21 

2.88 

8.61 

40.5 

21 

-  24 

2.17 

5.98 

33.5 
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the  BLOHARD  routine  for  both  the  approximating  beta  and  the  actual 
discrete  wind  speed  distributions  are  listed  in  Table  3.3-12  for  various 
costs  of  purchased  electricity.   There  is  little  difference  between  the 
optimum  WTGS  for  either  wind  speed  distribution.   For  low  electrical 
energy  costs,  an  optimum  WTGS,  which  yields  a  negative  annual  net 
savings,  is  found.   This  occurrence  indicates  a  local  optimum  for  which 
the  amount  of  money  lost  by  using  a  WTGS  is  minimized  although  the  true 
optimum  would  be  a  zero-size  WTGS. 

From  Fig.  3.3-17,  which  shows  the  WTGS  size  plotted  as  a  function 
of  electrical  energy  costs,  the  shape  of  the  trajectory  has  many 
characteristics  of  those  found  for  the  single  and  dual  season  idealized 
wind  speed  and  load  models.   For  very  high  costs  of  electricity,  the 
rated  power  of  the  optimum  WTGS  tends  toward  the  maximum  load  for  the 
year.   With  decreasing  electrical  costs,  rated  power  decreases  some- 
what but  rated  speed  increases  at  a  faster  rate.   Finally,  as  electrical 
costs  approach  breakeven,  the  rated  power  begins  to  decrease  at  a 
faster  rate  than  that  at  which  the  rated  speed  increases.   Although  this 
realistic  case  does  not  exhibit  the  dramatic  step-like  changes  obtained 
for  the  model  wind  speed  and  load  profiles,  the  overall  behavior  is 
still  quite  similar  to  these  model  cases. 

Breakeven  costs  can  be  determined  from  a  plot  of  annual  net  savings 
versus  cost  of  electricity  as  shown  in  Fig.  3.3-18.   For  this  farming 
enterprise,  the  cost  of  electricity  must  be  about  6.8c/kWh  in  order 
to  offset  the  installation  cost  by  the  savings  in  purchased  electricity. 


131 

>, 

1           O 

u-i     C 

<H     0) 

en 

CO' 

rn 

CO 

r-^ 

CO 

en 

QJ    -H 

en    u 

en 

CO 

r^ 

r» 

ON 

CN 

CO 

s* 

en 

en 

CN 

H 

en 

en 

&-S    "4-1 

U-i 

0! 

4J 

0)    tc 

Z    o£ 

C    /-. 

o 

en 

On 

■* 

m 

CN 

— i 

en 

i-}  «rt  <* 

m 

r- 

t£ 

O 

m 

CN 

CO 

vD 

TO     >  w 

OTl 

vD 

cn 

in 

CN 

3     TO 

CN 

H 

!-H 

1 

1 

C     CD 

4J 

c 

CD 

< 

CD 

J= 

<f 

<r 

<r 

en 

en 

<r 

<T 

r-> 

3 

O 

o 

O 

o 

o 

o 

o 

O 

tC       •   <-v 

H 

i— i 

U 

3      4-1     ,£ 

CD 

H     U  3 

X 

x 

X 

X 

X 

X 

X 

X 

4-1 

C    (Di 

e 

l_,    ,— i  s^- 

vO 

r^ 

to 

o 

CC 

vO 

c 

r- 

•H 

3     QJ 

in 

ON 

en 

cc 

r- 

C 

-T 

CO 

£ 

W 

CN 

On 

CN 

CN 

fN 

CO 

ra 

m 

C 

IT| 

m 

m 

« 

m 

in 

fD 

O 

c 

o 

o 

C 

o 

o 

o 

tf 

•  -a  *-» 

H 

H 

H 

i-i 

W    0)  JS 

0  -3   2 

X 

X 

X 

X 

X 

X 

X 

X 

0 

0)     0)   ^ 

IM 

T— I        QJ     ^ 

on 

r~- 

<r 

CO 

un 

CO 

fa    c 

CO 

VO 

o 

<r 

o 

CN 

CO 

o 

cc 

ON 

H 

ON 

rH 

■H 

(J 

en 

■H 

•a 

-<r 

<r 

O" 

-d- 

<r 

*ff 

<t 

<r 

M 

a 

o 

o 

o 

c 

o 

q 

o 

o 

CD 

cu  tj  r-s 

rH 

H 

<U 

K  2 

cj 

X 

X 

X 

X 

X 

X 

X 

X 

fD 

•  p  Jjd 

N 

4J             "-^ 

m 

CO 

r*. 

ON 

en 

in 

TO 

O    i£ 

on 

on 

CN 

r» 

CN 

cc 

QJ 

U 

fa 

«c 

U"1 

-3- 

<r 

CN 

VD 

m 

CN 

4J 

3 

a 

4J 

ag 

•              4-) 

<? 

\c 

<r 

vO 

CN 

CN 

c 

O 

C       U  0 

0)    >     (= 

U"l 

tO 

cc 

ON 

fcO 

CO 

rH 

*o 

u       ^ 

H 

H 

CN 

e 

TO 

4) 

a 

N 

3          /-. 

(%. 

1/1 

m 

en 

PI 

■j? 

1** 

en 

•H 

£      ^3 

to 

■rt  B-  Jjj 

vD 

fcfl 

•$ 

CO 

SO 

m 

*» 

o 

4-i             •—• 

CN 

CN 

CN 

CN 

H 

CN 

CN 

r-\ 

w    g 

Q. 

cj>    U 

© 

H    TO 

3  fa 

0i 

e  a 

4J 

3     3 

M-l     -r-|   ^ 

B  -= 

°      U    T= 

-H   5 

■H     DC 

CN 

c 

CO 

Ps 

vC 

o 

CO 

vO 

JJ    M 

4J     }-.  Jjj 

rH 

o 

C 

O 

O 

c 

a  o 

«     4-1  *-*, 

O   W 

0     V  <J> 

Ed 

C 

© 

o 

o 

C 

o 

o 

c 

CN 

1 

1 

en 

p 

en 

■0    c 

14-1     0)     0 

aj 

o    ai  «h 

(1) 

4-1 

rH 

a  4-i 

,0 

11   CO     3 

TO 

CJ 

u 

TO 

C         J2 

4-1 

H 

Ty 

Wind 

Distri 

CD 
P5 

0 

CO 

^ 

2  h 

o    J 

w  1 

^ 

o 

o 

"1 

o    y 

o 

^> 

cc 

O 

d 

-        / 

o      / 

O  ^ 

0) 

d 

-o 

a 

« 
n 

01 

g 

1-1 

a 
a 

e 

•H 

1 

•rl 
U 

i-l 
M 

4J 

a) 

V 

4J 

en  *h 
*J    u 

01    v 

(A     4-» 

HI     U 
(J     01 
G.H 
O     0) 
H 

iM 
r-l     0 
0 

&    u 

E   w 
^  o 

en    u 

w    d 

CO    o 

- 

to 

U-t 

c   > 

m    1 

M 

0 

U 

c    I 

g 

10 

0>     tc 

a  a) 

d  P 

o 

a  u 

! 

3  XI 

1 

i 

(suou^)    paads   psi^ 


133 

1800 

1600 

1400 

- 

^ 

</> 

1200 

— 

in 

oo 

■H 

> 

1000 

CO 

to 

JJ 

Cv 

2 

800 

r^ 

ce 

s 

600 

< 

400 

200 

0 

200 

1            1           1            1 

6 

7          8          9         10 

9 

Cost  of  Electricity  ($/kWh  x  10") 

Fig- 

3.3-18 

Annual  net  savings  versus  cost  of  electricity 
for  optimum  WTGS  -  Case  14. 
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This  value  is  subject  to  change  with  changes  in  the  interest  rate  and  the 
lifespan  of  the  WTGS.   As  mentioned  earlier,  an  interest  rate  of  10% 
per  year  and  a  lifetime  of  20  years  were  assumed.   These  values  are 
typical  for  today's  economy  and  used  in  order  to  get  an  indication  of  the 
feasibility  of  the  WTGS.   For  more  accurate  results,  more  precise  values 
of  interest  rates  and  useful  life  costs  should  be  used,  as  well  as  more 
accurate  cost  models  for  the  type  of  WTGS  to  be  installed. 

For  this  realistic  example  and  with  the  stated  economic  assumptions, 
the  breakeven  cost  of  electricity  is  well  above  the  current  price  charged 
for  purchasing  electricity.   However,  conclusions  as  to  the  feasibility  of 
a  WTGS  are  dependent  upon  the  technical  and  economic  assumptions  made. 
A  lower  cost  model  for  the  WTGS  will  definitely  alter  the  results  towards 
making  the  WTGS  economically  more  attractive.   Similar  results  would 
be  obtained  if  low  interest  loans  become  a  reality  (e.g.,  through 
government  subsidies).   Decisions  concerning  these  economic  parameters 
and  questions  must  be  made  before  any  definite  statement  is  made  as  to 
the  economic  feasibility  of  the  WTGS. 

Although  the  results  for  this  realistic  case  show  the  WTGS  to 
be  currently  economically  infeasible,  a  more  important  result  is  the 
observation  that  wind  speed  and  load  data  must  be  supplied  on  a  multiple 
season  basis  in  order  to  get  an  accurate  value  for  the  optimum  size 
WTGS.   As  less  data  are  supplied,  the  more  inaccurate  is  the  optimization 
procedure  and  the  calculated  optimization  results.   The  interaction  of 
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one  season  with  another  is  as  important  as  the  interaction  of  diurnal 
variations.   Therefore,  for  the  optimization  methodology  to  yield 
accurate  results,  the  input  data  should  be  as  complete  as  possible  if 
the  methodology  described  in  this  chapter  is  to  be  used  for  analysis 
of  actual  wind  energy  applications. 
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4.   CONCLUSIONS  AND  RECOMMENDATIONS 

In  this  research,  two  important  aspects  of  the  extraction  of  energy 
from  the  wind  were  analyzed.   First,  as  an  aid  to  help  predict  the  power 
available  from  the  wind,  analytical  wind  speed  distribution  models  were 
examined  along  with  methods  to  estimate  the  model  parameters.   Three 
techniques  were  presented  for  the  estimation  of  the  parameters  of  a 
modeling  Weibull  distribution,  which  is  the  most  commonly  used  analytical 
representation  of  observed  wind  speed  data.   In  addition,  the  beta  distri- 
bution was  introduced  as  an  alternative  wind  speed  distribution  model. 
Two  goodness  of  fit  tests  were  performed  on  each  analytical  distribution 
to  test  the  appropriateness  of  each  model  in  describing  observed  wind 
speed  distributions.   Second,  a  methodology  was  described,  whereby,  given 
wind  speed  distribution  models  for  a  particular  location  and  power  demand 
data  for  a  particular  energy-consuming  enterprise,  the  economically  optimal 
size  WTGS  could  be  determined  such  that  the  net  annual  economic  savings 
realized  with  the  WTGS  are  maximized.   It  was  assumed  that  the  WTGS  produced 
utility-compatible  electric  power,  and  furthermore  was  connected  into  the 
utility  grid  so  that  when  the  WTGS  could  not  generate  all  the  required 
power,  the  deficit  could  be  purchased  from  the  power  grid. 

It  has  been  shown  that  for  parameter  estimation  of  the  Weibull 
distribution,  the  matching-moments  technique  yielded  a  Weibull  distribution 
that  represents  the  wind  speed  data  much  more  closely  than  least  squares 
fitting  techniques.   From  the  results  of  the  analyses  of  544  observed 
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2 
wind  speed  distributions,  the  x   goodness  of  fit  test  indicated  little 

as  to  the  accuracy  of  either  of  the  parameter  estimation  techniques 
(matching-moments  or  least  square)  in  modeling  most  of  the  observed 
wind  speed  distributions.   However,  the  large  x  values  obtained  for 
approximating  Weibull  distributions  fit  by  the  matching-moments  method 
were  found  to  be  almost  totally  a  result  of  poor  fits  at  very  low 
speeds,  i.e.,  below  four  knots.   At  intermediate  and  high  wind  speeds 
(the  regions  most  important  for  a  WTGS),  the  Weibull  matching-moments 
fits  generally  produced  excellent  approximations  of  the  observed  distri- 
butions.  The  beta  distribution,  whose  parameters  were  also  estimated 
by  a  matching-moments  technique,  exhibited  a  similar  result  -  poor  fits 
at  very  low  speeds  but  excellent  fits  at  intermediate  and  high  speeds. 

The  second  procedure  developed  to  test  the  fitted  Weibull  and  beta 
distributions,  the  power  ratio  test,  confirmed  the  observation  that  the 
matching-moments  parameter  estimation  technique  represented  the  wind  speeds 
greater  than  four  knots  very  accurately.   With  this  test,  the  ratio  was 
calculated  of  the  power  obtained  from  a  given  size  WTGS  when  the  analytical 
(fit)  wind  speed  distribution  was  used  to  the  power  generated  by  the  same 
WTGS  when  the  discrete  (observed)  wind  speed  distribution  was  used. 
Application  of  this  test  to  544  observed  wind  speed  distributions  showed 
that  the  matching-moments  technique  yielded  Weibull  and  beta  distributions 
with  power  ratios  that  were  very  close  to  the  ideal  value  of  unity.   The 
least  squares  techniques,  however,  produced  Weibull  distributions  that, 
when  the  power  ratio  test  was  performed,  yielded  discrepanies  of  as  much 
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as  70%  from  unity.   Consequently,  for  both  the  Weibull  and  beta 
distributions,  the  matching-moments  technique  of  parameter  estimation 
provided  distributions  that  accurately  model  observed  wind  speed  distri- 
butions. 

It  was  also  seen  that  because  only  two  wind  speed  statistics,  i.e., 
the  mean  wind  speed  and  variance  of  wind  speeds  about  this  mean  are  needed 
to  calculate  the  parameters  of  the  Weibull  or  beta  distribution  using  the 
matching-moments  technique,  the  need  for  detailed  historical  wind  speed 
information  is  eliminated.   Since  detailed  meteorological  wind  speed 
distribution  data  are  not  readily  available  for  most  locations,  much 
simpler  and  less  time-consuming  measurements  or  analysis  of  meteorological 
data  tapes  to  obtain  the  mean  wind  speed  and  the  variance  of  wind  speeds 
need  be  performed.   Although  both  Weibull  and  beta  distributions  give 
accurate  fits  of  the  wind  speed  distributions  when  the  parameters   are 
estimated  by  the  matching-moment  technique,  the  beta  distribution  is 
particularly  attractive  since  its  two  parameters  can  be  calculated  directly 
from  the  mean  and  variance  of  the  wind  speed.   Calculation  of  the  Weibull 
parameters,  on  the  other  hand,  requires  the  numerical  solution  of  trans- 
cendental equations . 

From  the  sensitivity  studies  performed  on  the  optimization  methodology 
in  part  two  of  this  work,  it  was  seen,  that  in  addition  to  supplying  de- 
tailed wind  speed  distributions  and  load  demands  to  characterize  accurately 
diurnal  variations  in  wind  speeds  and  load  requirements,  seasonal  variations 
need  to  be  represented  also.   If  only  a  single  season  is  used,  the  optimal 
WTGS  size  as  a  function  of  cost  of  purchased  electricity  does  not  vary 
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smoothly  as  expected,  but  exhibits  very  rapid  changes  in  size.   However, 
as  more  seasons  are  added  to  the  analysis,  the  trajectories  of  the 
optimal  WTGS  size  as  a  function  of  cost  of  electricity  vary  in  a  much 
smoother  fashion.   Besides  needing  detailed  diurnal  wind  speed  distri- 
bution and  load  demand  requirements,  seasonal  variations  in  these  two 
inputs  must  also  be  characterized.   However,  single  season  sensitivity 
analyses  are  useful  in  identifying  the  general  trends  of  the  optimization 
methodology. 

Of  all  the  sensitivity  studies  analyzed,  the  parameter  that  had  the 
largest  effect  on  the  optimum  WTGS  was,  as  expected,  the  mean  speed  at 
the  WTGS  site.   For  even  the  single  season  case,  locations  that  have  mean 
speeds  of  20  knots  were  capable  of  producing  electricity  on  a  competitive 
basis  with  utility  supplied  energy.   Another  parameter  which  affects 
significantly  the  size  of  the  optimum  WTGS  was  the  size  of  the  load 
served  by  the  WTGS.   For  large  loads,  more  money  can  be  saved  from  the 
installation  of  an  optimum  WTGS,  and  consequently,  power  can  be  produced 
more  cheaply.   This  preference  for  larger  loads  is  a  direct  result  of  the 
lower  unit  capacity  costs  for  larger  UTGSs  inherent  in  the  WTGS  cost 
model.   Although  the  unit  capacity  cost  ($/kWh)  generally  decreases  as 
the  rated  power  of  a  WTGS  increases,  a  power  level  is  eventually  reached 
(a.  40-60  kW)  above  which  this  observation  is  no  longer  true  and  the 
cost  ($/kW)  actually  increases.   This  is  caused  by  the  fact  that  for 
large  rated  powers,  the  complexity  and  sophistication  of  the  WTGS  must 
increase  compared  to  smaller  machines  in  order  to  handle  the  large  amount 
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of  power.   Consequently,  this  transition  power  must  be  clearly  noted  so 
that  a  different  cost  model  can  be  assumed  to  account  for  this  added 
system  complexity. 

Another  parameter  that  was  shown  to  have  only  a  slight  effect  on 
the  optimum  WTGS  was  the  relative  magnitude  of  the  peaks  in  the  load 
demand.   As  long  as  the  load  peaks  coincide  or  are  inphase  with  the 
daily  wind  speed  peaks,  only  marginally  more  net  savings  and  lower 
breakeven  costs  are  achieved  if  the  load  varies  smoothly  about  an 
average  wind  speed  than  a  load  with  large  or  flashy  variations  about 
the  same  average  wind  speed.   However,  a  constant  load  demand  throughout 
the  day  has  lower  net  savings  and  a  higher  breakeven  cost  than  the 
smoothly  varying  load  with  the  same  average  load  demand  as  the  constant 
load.   Hence,  load  peaks  whether  large  or  small,   which  are  inphase  with 
the  wind  speed  peaks  tend  to  favor  the  development  of  wind  power  by 
lowering  the  breakeven  value. 

The  relative  magnitude  of  diurnal  wind  speed  variations  were  found 
to  affect  the  size  of  the  optimum  WTGS  to  a  lesser  degree.   The  optimum 
WTGS  produces  a  greater  net  savings  and  lower  breakeven  costs  for  a 
wind  with  large  diurnal  variations  than  a  wind  with  small  diurnal 
variations  about  the  same  daily  mean  speed.   As  was  true  in  the  load 
variation  cases,  as  long  as  the  peaks  in  the  wind  speed  and  load  variations 
are  inphase,  the  magnitude  of  either  do  not  affect  the  optimum  WTGS 
significantly. 

The  variance  of  wind  speeds  about  a  mean  wind  speed  exhibited  a 
detrimental  and  also  an  advantageous  effect  depending  upon  the  cost  of 
purchased  electricity.   Wind  speed  distributions  with  small  variances, 
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i.e.,   highly  peaked  around  a  mean  speed,  had  large  net  savings  for  high 
costs  of  purchased  electricity,  but  the  net  savings  decreased  sharply 
with  decreasing  purchased  energy  costs.   The  breakeven  value  was 
higher  for  wind  speed  distributions  with  a  small  variance  than  the 
distributions  with  a  large  variance.   However,  the  self-sufficiency 
of  the  optimal  WTGS  was  enhanced  when  the  location  site  had  wind 
speed  distribution  with  a  small  variance.   Consequently,  if  a  greater 
self-sufficiency  is  required,  sites  with  almost  constant  speeds  should  be 
chosen,  whereas,  if  net  savings  is  to  be  maximized,  sites  which  have 
wind  speed  distributions  with  large  variances  are  more  attractive. 

Finally,  the  investigation  of  the  effect  of  credit  for  surplus 
electricity  on  the  optimum  WTGS  showed  the  optimum  rated  power  was 
most  affected.   The  WTGS  rated  power  increased  greatly  compared  to 
similar  cases  which  received  no  credit  for  surplus  power.   This  increase 
is  reasonable  since  as  more  surplus  power  is  generated ,  more  money  is 
made.   For  the  small  credits  which  are  likely  to  be  given  (if  any  credit 
is  given  at  all),  no  significant  effect  on  annual  net  savings,  self- 
sufficiency,  or  breakeven  value  is  realized.   However,  as  noted  earlier, 
if  the  WTGS  rated  power  increases  beyond  a  certain  critical  value,  a 
higher  cost  model  needs  to  be  applied.   Hence,  this  greater  WTGS  cost 
may  overshadow  any  profit  from  selling  surplus  power. 

Many  areas  for  further  study  have  emerged  during  the  course  of  this 
work.   For  instance,  when  modeling  observed  wind  speed  distributions,  the 
effect  of  bin  size  (speed  subintervals)  used  in  the  histogram  of  observed 
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wind  speed  distributions  should  be  investigated.   This  study  used 
rather  large  speed  subintervals  when  characterizing  wind  speed  data. 
To  reduce  the  bin  size,  analysis  of  the  meteorological  data  tapes  would 
be  required.   Furthermore,  the  bin  size  can  be  reduced  only  so  far  before 
statistical  fluctuations  (caused  by  the  paucity  of  wind  speed  data)  will 
begin  to  mask  actual  detail  in  the  wind  speed  distribution.   The  whole 
area  of  extracting  the  most  information  about  the  wind  speed  distribution 
from  finite  amounts  of  observed  data  is  a  very  important  one  if  accurate 
analysis  of  wind  energy  is  desired. 

Although  only  two  analytical  distributions,  the  Weibull  and  beta, 
were  examined  in  this  study,  many  more  can  be  tried.   For  example,  the 
beta-prime  distribution  could  be  used  as  a  fitting  distribution.   The 
beta  distribution  used  in  this  study  is  defined  over  a  finite  interval, 
while  the  beta-prime  distribution  is  defined  over  the  entire  positive 
axis  [29].   Other  distributions  worthy  of  further  study  are  the  gamma 
distribution  [9]  and  the  single  parameter  Rayleigh  distribution  [7,8,9]. 

The  optimization  methodology  investigated  in  this  work  considered 
only  several  basic  features  of  WTGS  to  obtain  an  optimum  match  between 
wind  energy  production  and  demand  load.   Extensions  to  this  model  could 
be  made  to  incorporate  many  additional  features.   For  example,  a  logical 
extension  would  be  to  assume  a  distribution  function  about  a  given  mean 
value  which  would  give  the  variation  of  the  demand  load  rather  than 
assuming  a  constant  or  mean  load  during  every  time  interval.   Such  a 
change  would  require  only  slight  modifications  in  the  optimization 
methodology. 
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An  additional  sensitivity  study  would  involve  the  investigation  of 
load  management.   In  this  study,  the  load  demand  peaks  were  assumed 
a  priori,    and  no  attempt  was  made  to  manage  the  load,  within  given  con- 
straints, to  improve  further  the  savings  afforded  by  the  optimum  WTGS. 
Such  load  management  as  part  of  the  optimization  procedure  would  require 
the  inclusion  of  linear  programming  techniques  into  the  present 
methodology. 

Because  the  availability  of  low  cost  loans  would  affect  the 
feasibility  of  a  WTGS,  the  exact  nature  of  such  loans  is  yet  another 
parameter  which  could  be  incorporated  easily  into  the  objective  function 
of  the  present  methodology.   Similarly,  much  more  detailed  economic 
models  for  WTGS  costs,  amortization  incentives,  tax  credits,  surplus 
WTGS  energy  credits,  etc.  should  be  investigated.   Such  investigation 
would  require  only  minor  modifications  of  the  present  methodology  which 
has  been  written  in  a  highly  modularized  fashion  to  facilitate  alterations 
in  components  of  the  methodology.   As  previously  noted,  the  high  wind 
speeds  favor  the  economic  viability  of  a  WTGS.   Because  wind  speeds 
generally  increase  with  height  above  ground  [3],  a  WTGS  at  a  greater 
height  could  potentially  produce  larger  amounts  of  power.   However, 
as  WTGS  height  increases,  the  tower  cost  also  increases.   Consequently, 
with  a  cost  model  for  tower  costs  and  a  model  for  increasing  wind 
speed  with  height,  the  optimization  methodology  could  be  easily  modified 
to  include  WTGS  height  as  a  variable  to  be  optimized  in  the  selection 
of  the  economically  optimal  WTGS. 

This  study  assumed  no  storage  of  excess  power.   Conceivably  the 
inclusion  of  battery  storage  could  improve  the  economic  attractiveness 
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of  wind  systems,  and  the  effects  of  battery  storage  should  be  investigated. 
Because  the  optimum  amount  of  battery  storage  capacity  can  be  expected  to 
be  very  sensitive  to  the  amount  of  surplus  energy  generated,  much  more 
detailed  wind  speed  information  is  needed  so  that  an  accurate  estimate  of 
optimal  capacity  can  be  made.   This  would  require  either  model  wind 
speed  data  for  every  day  in  the  year  or  correlations  which  describe 
how  wind  speeds  depend  on  earlier  wind  speeds. 

Finally,  this  same  optimization  methodology  used  in  this  work 
can  be  applied  directly  to  solar  energy  studies  of  electrical  energy 
generation.   The  wind  speed  distributions  would  be  replaced  by  solar 
insolation  distributions,  i.e.,  the  probability  a  particular  amount  of 
radiant  energy  will  strike  a  surface  in  a  given  time  interval.   The  WTGS 
response  function  can  be  replaced  by  the  solar  cell  response  to  the  radiant 
energy  that  strikes  it.   If  solar  collectors  are  used  to  produce  thermal 
energy,  then  storage  capacity  would  have  to  be  included  in  the  optimization 
procedure  (analogous  to  battery  storage  for  the  wind  energy  problem) . 
Hence,  the  optimization  methodology  used  in  this  study  can  be  expected 
to  be  applicable  to  the  selection  of  optimal  solar  energy  systems  which 
maximize  the  system  parameters  so  as  to  produce  the  maximum  economic 
savings. 

In  conclusion,  this  research  was  two-fold.   The  first  part  studied 
the  modeling  of  observed  wind  speed  distributions  by  analytical  distri- 
butions.  The  second  part  applied  the  wind  speed  models  to  match  the 
available  wind  power  from  a  WTGS  with  the  load  demand  requirements  to 
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compute  the  size  of  an  optimal  WTGS,  i.e.,  a  WTGS  that  saves  the  user 
the  most  amount  of  money  by  replacing  normally  purchased  power  with  that 
generated  by  a  WTGS.   This  type  of  detailed  investigation  is  needed  in 
order  to  determine  the  impact  wind  generated  power  can  have  on  future 
energy  policies. 
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APPENDIX  A 
The  Program  CURVEFIT 

The  code  CURVEFIT  calculates  the  parameters  of  analytical  wind 
model  distributions  (the  three  Weibull  distributions  and  the  beta 

distribution)  by  fitting  these  distributions  to  observed  discrete 

2 
wind  speed  distributions.   In  addition,  the  x  and  power  ratio 

statistics  are  generated  as  a  measure  of  the  goodness  of  fit  of  the 

analytical  distributions.   The  program  input  requires  the  roots  and 

weights  of  any  even  order  Gauss-Legendre  quadrature  desired.   For  the 

computation  of  the  power  ratio,  the  parameters  of  the  WTGS  (i.e.,  rated 

power,  rated  speed,  cut-in  speed)  are  needed.   Finally,  the  boundaries 

of  the  wind  speed  subintervals  used  in  the  discrete  distribution  are 

required  along  with  the  frequency  of  wind  speed  observations  in  each 

subinterval  and  the  total  number  of  wind  speed  observations  made. 

Detailed  input  requirements  are  described  by  comment  cards  in  the 

FORTRAN  IV  program  listing  which  is  included  in  this  Appendix. 

Estimation  of  the  parameters  of  the  fitting  distributions  is 

performed  according  to  the  methods  described  in  subsections  2.3.1  and 

2.4.1.   The  numerical  procedure  used  to  solve  Eq.  (2.3-18)  for  the 

parameters  of  the  Weibull  distribution  using  the  matching-moments 

method  is  the  subroutine  RTMI  [30],   The  subroutines  WBLFIT  and  F 

compute  the  initial  points  for  the  Mueller's  iteration  technique  used 

in  RTMI. 
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The  calculations  required  for  the  x  and  power  ratio  tests 
described  in  Section  2.5  are  performed  by  the  subroutine  CH1SQ  and 
POWER,  respectively.   Any  Weibull  or  beta  distribution  values  needed 
by  CHISQ  and  POWER  are  computed  by  the  subprograms  WB  or  FI,  respectively. 
The  necessary  integrals  are  computed  by  the  subroutine  GAUSS  which  uses 
a  specified  even  order  Gauss-Legendre  quadrature . 

The  output  of  CURVEFIT  is  either  a  table  of  model  distribution 
parameter  values  or  a  table  of  goodness  of  fit  statistics  (similar  to 
Tables  2.6-1  and  2.6-2)  depending  upon  the  values  of  program  option 
parameters.   The  program  is  written  in  FORTRAN  IV  for  the  Kansas  State 
university  ITEL  Advanced  System  5  (which  is  equivalent  operationally  to 
an  IBM  370/158  system).   Liberal  use  of  comment  cards  and  variable  names 
with  high  mneomic  content  assist  the  user  in  deciphering  the  program 
logic. 
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FORTRAN    IV    G    LEVEL      21 


DATE    «    7E135 


***** 


CURVEFIT    **»* 


****** 


********** 


THIS    PROGRAM    CCHPUTES    THE    PARAMETERS    OF   THE    WEIBULL    CI STRI BUTICNSI LEAST 
SCUARES,    WEIGHTED    LEAST    SQUARES,    AND    HATCHING-MGHENT5    TECHNICJES)     ANC   THE 
BETA    DISTRIBUTICNIMATCHING-HCHENTS   TECHNIQUE].      CHI-SQUARE    AND    POWER   RATIO 
STATISTICS    ARE    ALSO   COMPUTED. 

INPUT   DATA! 

CARD    1  FORMAT    (12) 

NHALF    =    THE    HALF    VALUE    OF    THE    EVEN    ORCER    GAUSS-LEGENDRE    QUADRATURE 
USED    TO    EVALUATE    THE    NECESSARY    INTEGRALS 

C*  CARD    2  FORMAT    (4G20.7) 

C*  ROOTII)    "    QUADRATURE    ORD1NATES (ONLY    POSITIVE    VALUES) 

C*  (MAY    BE    MANY    CAROS) 

C» 

C*  CARD  3  FORMAT     (4G20.7) 

C*  WEIGHT! I  I     =    QUADRATURE    WEIGHTS 

C*  (MAY    BE    MANY   CARDS) 

C* 

C*  CARD    4  FORMAT    (415) 

C*  INN  *    TOTAL    NUMBER    OF    POSSIBLE    SPEED    SUBI NTERVALS (USUALLY    ID 

C*  IPR    *    RATED    POWER    OF    WIND   TURBINE(KH) 

C*  IVR   =    RATED   SPEED   OF    WIND   TURBINE(MPH) 

C*  IVC   ■    CUT-IN    5PEEO    OF    WIND    TURBINEIMPH] 


OOOl 

0002 
0003 
0004 
0005 
OOOt 

0007 

0008 


C* 
C» 
C* 
C* 

c« 
c* 


MM!  5     rO^'AT[12FS.2l 

IVINT(I)  -  ENDPCINTS  OF  WIND  SPEED  SUE  INTERVALSIKHOTS) 

CARD  6     FCRMATI20A4) 

TITLE  -  TITLE  CARD  FOR  WIND  SPEED  DATA  SET  ANALYZEO 
(HAY  BE  MANY  CAROS) 

CARD  7     FORMAT! 12, 13,1014,16,151 

MONTH  *  MONTH  FROM  WHICH  KIND  SPEED  DATA  IS  OBTAINED 

NTIME  =  DAILY  TIME  PERIOD  FROM  WHICH  WIND  SPEED  DATA  IS  OBTAINED 

IFREQ(I)  =  FREOUENCYtX  1000)  OF  OBSERVATIONS  If.  I-TH  SPEED  SUBINTERVAL* 

(BEGINNING  WITH  FREQUENCY  IN  2ND  SPEED  SUBINTERVAL) 
NUMKY  *  SPACE  FOR  DATA  I  CENT  IF ICAT ICN  PURPCSEStCAN  ALSO  BE  USED  TO 

ADD  TWO  MORE  SPEED  SUB  INTERVALS) 
NOBS  =  TOTAL  NUMBER  OF  WIND  SPEED  DATA  OBSERVATIONS 
(HAY  BE  MANY  CARDS) 


C» 

C* 

c» 

C» 

c* 

c* 

C*  WRITTEN  BY  L.  A.  POCH,  KANSAS  STATE  UNIVERSITY,  DECEMBER  1977. 
C» 


*************** 


P.EAL*4    FREC(41l,TITLE(40l,VINT(4i),Y(41],X(41),'V(4I),FRECN(41)t 
1    F(41),RCGT(20] , WEI GHT ( 20  I , FREYl 41) 
REAL* 4    MMEAN,IVINT<41) ,K,MEAN 
iNTEGEP*4    IFREC141) , DF1 ,DF2 ,DF3 ,DF4  '     . 
CCMMCN/LINK1/K,C,AA,BB,VMAX,  VRATED.FCTR 
COMMON/ LI  NK3 /NHALF, ROOT, WEI GHT 
EXTERNAL    W6 ,F I , V3F I , V3WB 
«*»    READ    IN   GAUSS-LEGENDRE    QUADRATURE    ORDINATES    AND    WEIGHTS 
READ(5,112INHALF 
113    FCRMATII2)  '     ' 
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FORTRAN  IV  G  LEVEL   21 


MAIN 


DATE  "  78135 


23/28/52 


0  009 
0010 
0011 

0012 
0013 
0014 
0015 

0016 
0017 

0018 
0019 
0020 
0021 
0022 
0  023 
0  024 


0025 
0026 
0027 


0028 
0029 


0030 
0031 


0033 
003* 
0035 
0036 
0037 
0038 

0  039 
0040 
0041 
0042 


PECIFICATIONS 


READ (5, 12}IR00T(I),I=1,NHALFI 
READ t 5, 12 M WEIGHT  (II  , I*1,KHALFI 

12  FORMATI4G20.7I 
C***    READ    IN    NUMBER    OF    SPEED    SUSINTERVALS    AND    WIND    TLRBINE 

READ(5,2)INN,IPR,IVR,IVC 

2  F0RHAII4I51 
NI\*INN*1 
PRATED«=IPR 

C***    CCNVERT    MPH    TO   KNOTS 

VRATED=IVR/1.15 

VCUTIN-IVC/1.1S 
C**»    READ    IN    ENDPUNTS    OF    SPEED    SUSINTERVALS   AND    TITLE 

READ! 5t II  (  IVINT(I) , I*1,NN> 
1    FORMAT! 12F5. 21 

3  PEAD(5,100)    TITLE 

100  FORMATI20A4) 

WRITE  16, 110 1 (TITLE  I  11,1-1,331 
110    FCRMAT!'1'////,33A4I 
PRINT    1001 
C--*    REMOVE    COMMENT    TO    LIST    CHI-SQUARE    AND   POWER   RATIO   TABLE 
C1001    FORMAT! '0' ,1181 •-') J 
1001    FORMAT! 'C  ,113!  •-■)) 
PRINT    102 
102    FORMATIT52. 'WEIBULL    DISTRIBUTION    PARAMETERS  * ,T98, 'BETA    D1STRIBUTI0 
IN'/'    MONTH      TIME         MEAN    SPEED         SID.    OEV.         LST.    SOS.-UNWTD.         LST 
2.    SOS.-WTD.        HATCHING-MOMENTS  PARAMETERS'/TV, 'IHRSI         (UNO 

3TSI  (KNOTS! ' ,3I6X,'K  C    'J, EX, 'ALPHA  BETA') 

C«*»    REMOVE    COMMENT    TO   LIST    CHl-SCUARE    AND    POWER    RAT1C   TABLE 
C    102    FCRMATIT27, 'RESULTS    OF    CHi-SOUARED    TEST ' ,T8 1, 'RESULTS    OF    POWER    RAT 
C  1IC    TEST*'/T24, 'WEIBULL    DI STRI BUTICN* ,T79 , 'WEI  BULL    DISTRIBUTION'/ 

C  2     •    MCNTH      TIME         LST.    SOS.  LST.    SOS.  MATCHING- ', 7X, ' BETA' , 

C  3    9X,'LST.     SCS.         LST.    SCS.         M4TCMI NG-' .7X, ' BETA'/ 

C  4T9,  MHRS)  IUNWTC.) ' ,6X, '£ WTO. I ', 6X ,' MOMENTS' ,5X , 'DISTRIBUTION' , 

C  5    6X, • (UNWTC.)'  ,5X,' IWTD.)  ' , 4X , 'MOMENTS' ,5X,  'DISTRIBUTION »1 

PRINT    1001 
99    N=INN 
C»«»    READ    IN    FRECUENCY    OF    OBSERVATIONS 

READi5,101,END-98)HONTH,NTIME,(lFREC(I>  , 1-2,N|  ,NUMMY,    NOES 

101  FCRMATCI2,I3,10I4,I8,15I 

C-"    IF    MOPE    THAN    CNE    LOCATION    IS    TO    BE    ANALYZED    A    BLANK    CARD    FDR    DATA    CARD   7 
C  WILL    CAUSE    THE    PROGRAM    TO    INCREMENT    TO    THE   TITLE   CARD    OF    THE    NEXT    DATA    SET 

IFINOES    .NE.    0)    GO    TO    13 
C»**    REMOVE    COMMENT    TO   LIST    CHI-SQUARE    AND    POWER    RATIO    TABLE 
C  PRINT    997, IPR,IVR,IVC 

C    997    FORMAT! '-»      POWER   RATIO    COMPUTED    FOR    RATED    POWER    *    ',13, 
C  1    '    KW,    RATED    SPEED   =    ',12,'     MPH,    CUT-IN    SPEED   «    ',12,'    MPH'/ 

C  2    •    ',1181 '-"II 

GO   TO   3 

13  DO    172    I-l.NN 
172    VINTIII-IVINTUI 

DC   401    1=2, N 
401    FREQ(II«IFREO(I>/1000. 
SUM=0.0 
C-»*    COMPUTE    FREQUENCY    CF    OBSERVATIONS     IN    INITIAL    SPEED    SUBINTERVAL 
DO    54    l»2,N 
54    SUM»SUM*FREOIII 
FRE0.I1I-1.0-SUN 
SUM«SUH»FREO(l) 
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FORTRAN    IV    G   LEVEL       21 


HAIN 


DATE    -    7S135 


0043 
0044 
0045 
0046 
0047 
0  048 
0049 
0050 
0051 
0052 
0053 

0054 
0055 
0056 
0  057 
0058 
0059 

0060 
0061 
0  062 
0063 
0064 
0065 
0066 


0067 
0066 
0  069 
0070 
0071 
0072 
0073 
0074 
0075 
0076 
0077 
0078 

0079 
OOSO 

0081 
0082 
0083 
0064 
0085 
0086 

0087 
0088 


c».»    COMMUTE    DISCRETE    CUMULATIVE    DISTRIBUTION    FUNCTION 
SUK1-0.0 
00   41    I-l.N 
FRE0(1I»FRE0(1I/SUN 
SL'».l»SUfll»FRECm 
41    FCII-SUH1 

DC    130    1-1, N 

I F I F 1 1 1     .GE.    .999999JG0   TO    131 

130  CONTINUE 

131  N-I 
VMAX-VINTIN+11 
SUM2=0.0 

C...    COMPUTE    DISCRETE    PRDBAEILITY    DENSITY    FUNCTION 
DO    46    IK=1,N 
SUK2  =  SU.W.2*FREQ(IK) 
FIIK1-SUK2 

Y1IK1-VINT! IK<1I-VINTCIK) 
FRECM1K)  =  FREG(IK)/YCIK) 
46    V1!M=0.5>I  VINTUK+11+VINTI  IKH 
C««»    COMPUTE    MEAN,    VARIANCE,    ANC    STANDARD    DEVIATION    CF    NINO   DAT* 
HEAN-0.0 
VAR-0.0 
DC   30    1=1, N 
3C    HEAN=MEAN+VII1*FREQIII 

DO    31    I-l.N 
21     VAR-VAR    *     (V(II-NEAN)* 
STDEV-SOP.TIVAR) 
c».»    COMPUTE    WIND    TURBINE    GENERATOR    POWER    USING    DISCRETE    PROBABILITY    DENSITY 
C  FUNCTION 

HPOWER=0.0 

DC   300    1-1, N 

IFtVINTU+ll     .LT.    VCUTIN)    GO    TO   300 

IF(VINTII)     -GT.    VRATEDI    GO    TO    304 

VCUTIN  .AND,  VINT (1*1 1 
VINTII+11 
V I  NT [ 1*11 
VINT  11  +  1 1 


»2»FRE0(II 


.  LT. 


VCUTIN  .AND. 
VCUTIN  .AND. 
VCUTIN  .AND. 


.GT.  VRATEDI  GO  TO  302 

.LT.  VRATEDI  GO  TO  303 

-LT.  VRATEDI  GO  TO  305 

.GT.  VRATECI  GO  TO  333 


332 


IFIVINTI I  1 
IF(VINTU) 
I F I V 1  NT (I) 
IFtVINT! I) 
PRINT    332,1 

FCRMATl'OINTERVAL    DOES    NOT    FIT    ANY    CATEGORY ', 5>, ■ INTERVAL"" ,151 
GO    TO    300 

333    FFDWER-HPCWER-MC.25»(VR1TED»»4-VINTII>»»4)/VRATED»«3   ♦ 
1    VINTU+W    -    VRATEDI*FREO.N(II 
GO   TO    300 

302  HPCWEP-MPOWER    ♦    10.25*1 VRATED**4   -    VCUT IN**4] /VRATE0**3    + 
1    VINTU+1)    -    VF.ATED1-FRECMI] 

GO   TO    300 

303  hPDWER-HPOWER    »    0.25»FRECN(  1 1  <■(  VINTt  1*1 )  «»4   -    VCUTIN»»4! /VRATE0»»3 
GO    TO    300 

304  HPO;iER=HPOHER    «    FREONI 1  l»l  VINT1  I*1I-VINTI  II  ) 
GC   TO    300 

305  HPCHER-HPCHER    »    0.25»FRECM 1 1 •( VI  NT  1 1*1 )*»4    -    VINTII)»«4I/ 
1    VRATE0**3 

300    CONTINUE  ' 

HPCWER=PRATED»HPOWER 
0 

C***    WEIEULL    FITTING 
C»»»    LEAST    SOUARES    SOLUTION  FOR    THE    HEIBULL   FIT 

SUNX-0.0 
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0090 

SUHY=0.0 

0091 

SUMXY=0.0 

0092 

SUMX2=0.0 

0093 

N1=N 

009* 

DD    20    1=1, Nl 

0095 

YU]=ALGG(-ALOG( l.OOOl-FC 1) )J 

0096 

X(I.=ALGG(V(I>) 

0097 

SUMX2=SUMX2    ♦    X(I)**2 

0098 

SUMX=SUMX    +    Xt I) 

0  099 

SUMY=SUMY    +■    Y(II 

0100 

20 

SUMXY=SUMXY    +    X(1)*YIIJ 

0101 

K=(SUMXY    -    SUMX»SUMY/Nl I/.SUMX2   -    SUMX**2/Nll 

0102 

6=(SUMY    -    K*SUMX)/N1 

0103 

C=EXP(-B/K> 

0104 

UWTK=K 

0105 

UWTC=C 

0  136 

MMEAN=C*GAMMA( 1.0+1.0/K) 

0107 

VVAR*C=*2aGAMNA(  WQ  +  2.0/K)    -    MMEAN*»2 

0108 

FCTR=1.0 

C* 

.* 

CCMPUTE    CHl-SCUARE    AND    POWER    RATIO    STATISTICS 

0109 

CALL    CHISCIVINT,FREC,W&,NCBS.CHI2L.N,2,DFlf 

0110 

c 

c» 

CALL     POWER  1  PRAT  ED  ,HPOWS  R  ,  VRATED  ,  VCUTI  l\  ,  VMAX  ,V3hfi  ,  WB  ,WP1  iCAPFCll 

»* 

WEIBULL    FIT    WITH    A    LEAST    SQUARES    WElGhTEC    BY    OBSERVED   FREQUENCY 

0111 

SUMX-0.0 

0112 

SUMY* 0.0 

0113 

SUMXY-0.0 

0114 

SUMX2-0.G 

0115 

D3   iO    I-»li'i 

0116 

Y(I)=ALOG(-AL05(1.0001-F( I))) 

0117 

X(II*ALGG(VU)) 

011S 

SUMX2=SUf'.X2    +    FREGI  I)*XU)**2 

0119 

SUMX=SUMX    +    FREQ(n*XIIl 

0120 

SUHY=SUMY     +    FREQtI)*Y(l) 

0121 

6C 

SUMXY=SUMXY    +    FREQI I)*X(I )*Y( I) 

0122 

K=(SUMXY    -    SUMX«SUMY)/(SUMX2    -    SUMX**2) 

0123 

B=SUKY    -    K*SUHX 

0124 

C=EXP(-S/K» 

0125 

WTK*K 

0126 

WTC-C 

0127 

HM=AN=C* GAMMA, 1.0+1.0/K) 

out 

VVA«=C**2»GAHMAI 1.0+2. 0/K)    -    HMEAN**2 

0129 

FCTR=1.0 

c» 

** 

COMPUTE    CHI-SQUARE   AND    POWER    RATIC    STATISTICS 

0130 

CALL    CHlSC(VINT,FREQ,WB,NG3SrCHI22.N,2,DF2) 

0131 

CALL    POWER  (  PRATED,  HPOwER, VRAT  ED.  VCUTIN, VM4X.V3WB, W8, WP2 .CAPFC2) 

c 

c» 

*• 

CALCULATION    BY    MATCHING    MEAN   AND    VARIANCE    TO    KEISULL    DISTRIBUTION 

0132 

CALL    W3LFIT(HEAN,VAR,C,K,IER> 

0133 

IFUER.NE.01    GO    TO    35 

0134 

XHKK-K 

0135 

XHMC*C 

0  136 

MMEAN=C*GAMMAI 1.0+1.0/K) 

0137 

VVAR=C**2*GAMKAll.O*2.0/K)    -    NMEAN**2 

0138 

FCTR  =  1.0 

c» 

ft* 

COMPUTE    CHI-SQUARE    AND    POWER    RATIO    STATISTICS 

0139 

CALL    CHISQ(VlNT,FRE6iWB,NCBSrCH123,Ni2,0F3I 

0140 

CALL    POWEMPRATED,hPC.iER,VRATED,VCUTIN,Vf'AX,V3W5,WB#WP3»CAPFC3J 
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01*1  GO   TO    50 

0142  35    PRINT    36.IER 

01*3  36    FCRMATC     FIT    BY    HATCHING    MEAN    AND    VARIANCE    COULC    NOT    EE    ACHIEVED. 

1      RETURN    ERROR    CODE*' #131 

C 

C***     FIT    DATA    TO    A    BETA    DISTRIBUTION    BY    MATCHING    MEAN    AND    VARIANCE 
014*  50     VMAX*VINTIN+1) 

01*5  AA=<MEAN/VKAX1*(MEAN*(VKAX-MEAN)/VAR    -    l.C) 

Oi*6  BB=(VMAX-MEAN]»AA/MEAN. 

01*7  r*MEAN=AA*VKAX/(AA+BB) 

0148  VVAP=AA*fiB*VMAX**2/{ ( AA*BE) **2* ( AA*BS+1 .0 1 » 

01*9  FCTR=GAMMA(AA    +    BB 3  / ( VMAX*G AKMA ( A  A) »GAHKA tBBl ) 

C***    COMPUTE    CHI-SCUARE    AND    POWER    RATIO    STATISTICS 

0150  CALL    CHISCtVINT.FREQ.fi  ,NCBS  ,CH12* ,N ( 2, DF*» 

0151  CAlL    POWER  I  PRATED, HPOwER, VRATED, VCUT IN, VMAX ,V3FI .FI , HP*, CAPFC4I 

0152  NTIME=NTIME*1 

0153  NTIME1=NTI ME-3 

C*»»    REMOVE    COMMENT    TC   LIST    CHI-SOUARE    AND    POWER    RATIC    TABLE 

C  PRINT     160,MCNTH,NTIMEI .NT1ME.CF1 ,CHI2 1 , CF2, CHl 22, CF 3 , CHI  23. Of 4 

C  1     , CH I  24, CAP FC1, CAPFC 2, CAPFC 3. CAPFC* 

C    160    FORMAT  1 1 4, 3X.  1 2.  •-',I2.3X,4C' I  M2, •  > ', G9  .3 1 , 5X,G 10. 3 ,2X,G10. 3, IX, 

C  1    G10.3.4X.G10.3I 

0154  PR  TNT     160,MONTH,NTIMEl,NTIME,MEAN,STDEV,UhTK,UVTC,WTK,HTC,XHHK, 
1       XMMC,AA,BB 

0155  160    FCRMATU4,3X,I2.'-'  . 12, 4X , F7-3 ,6X , F7.3 , IX  ,3 13X , F7 .3 , IX, F7.3) ,3X, 

1    2F9.3) 

0156  GC   TC   99 

C**«    REMDVE    COMMENT    TO   LIST    CHI-SOUARE    AND    POWER    RATIC    TABLE 

C      98    PRINT    97,IPR,1VR,IVC 

r       97    FHSM4T(  >-*       POWPR    RATIO    CPMPIITED    FHR    RATFC    POWER    *     *.I3, 

C  1     ■    KM,    RATED    SPEED    =     '.12,'    MPHi    CUT-IN    SPEED    *    *,I2,'    MPH«  / 

C  3    '    '  ,118('-'l/'l*i 

0157  98    PRINT    97 

015B  97    FORMAT (•     ■ , 113 ( '-' 1 /'l1 I 

0159  STOP 

0160  END 
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0002 
0003 
0004 
0005 
0006 
0007 
0008 
0009 
0010 
0011 
0012 


FUNCTION    GICUAO    (A.BtFN) 
••    GAUSS-LEGENDRE    QUADRATURE    CF    FUNCTION    FN    CVER    INTERVAL    (A,B1 
«*    INTEGRAL    IS    SET    TC    1ERG    IF    LOWER    LIMIT    LARGER    ThAN    UPPER    LIMIT 

REAL**    ROOT120] ,WE1GKTI20> 

COPMCN/  LINK  It  NH  ALE  i  ROOT,  HEIGHT 

GLOUAD.0.0 

IF    IA.GE.8)    RETURN 

EA«0.5»(B-AI 

AB=0.5*IA*B1 

00    10    I-l.NHALF 
10    GLGUAD=GLQUAD+WEIGHT(I)»IFNtAe+8A»R00Tt  I)  >  +  FNUB-8A»RC0T(  11  I) 

GLOUAD=BA*GLO.UAO 

RETURN 

EHO 


0001  FUNCTION    F1(V> 

C«»«  SUBROUTINE  CALCULATES  VALUES  GF  EITHER  1 V/VRATEC)»*»3»EETA  OR  JUST  BETA 

C  DISTRIBUTION 

0002  COMHON/LINKl/K,C.AIJiEIJ,VMAXIJ,VRATEG.FCTR 

0003  FI  =  IV/VHAXIJ)»(AIJ-1.)»I1.-V/VHAXIJI»«1B1J-1.1 

0004  FI«FCTR»FI 

0005  RETURN 

nn(16  ENTRY   V3FKV1 

0007  V3FI=(V/VRATE0I»«3»IV/VHAXIJI»«IAIJ-1.I»(1.-V/VHAXIJ)»»(B1J-1.) 

0008  V3FI*FCTR*V3FI 

0009  RETURN 

0010  ENO 


0002 
0003 

0004 
0005 
0006 


SUBROUTINE    POWER  ( PRAT  EO.HFOWER,  VRATED,  VCLTIN,  VKAX.V3FC,  FC  ,  GEN, 
1    CAPFCI 

SUBROUTINE    COMPUTES    GENERATED    WIND    TUREINE    PCWER  FROK    ANALYTICAL 
DISTRIBUTION.       ALSO    COMPUTES    POWER    RATIO 

EXTERNAL    V3FC.FC 

GEN=PRATED*(GLQUADIVCUTIN,AKIN11VRATED,VKAXI iV2FCI* 
1    GLQUA01 VRATED, VMAX,FC) ) 

CAPFC=GEN/HPOW£R 

RETURN 

ENO 


159 


FORTRAN    IV    G    LEVEL      21  DATE    =    78135  23/28/52 

C 

0001  SUBROUTINE    CHISCI Z. FREQ, FT, J, CHI, N, ESTPAR ,DF) 
C«*»     SUBROUTINE    COMPUTES    CHI-SCUARE    VALUE 

C  ESTPAR    •    NO.    OF    PARAMETERS    OF    A    DISTRIBUTION    ESTIMATED   FROM    THE    DATA. 

C  DF    =    DEGREES    CF    FREEDOM    OF    CH-SCUARE    DISTRIBUTION 

0002  REAL**    Z(41)  .FREQKll  ,FCT(4l)  .FREY14U 

0003  INTEGER**    DF, ESTPAR 
000*  EXTERNAL    FT 

0005  CHl'0.0 

0006  CO   400    1=1, N 

0007  FREY(I)=FRECU) 

0C08  403    FCT( I )=GLQUAD(2( I) ,2(1*1) ,FT) 

0009  1N0EX1-0 

0010  INDEX2*Q 

0011  DC    401    II*ltN 

0012  JJ=N*1-II 

0013  XNUM=J*FCT( JJ) 

0014  IFIXNUM    .GE.    1.01    GO   TO    401 

0015  IFUJ    .EO.     1>    GO    TO    402 

0016  INDEX2=INDEX2*1 

0017  FCT(  JJ-II=FCT(JJ-1)+FCTUJI 

3  018  FFEY(  JJ-l)=FF.£Y(JJ-l)*FREYt  JJ) 

0019  NNN=N-JNDEXi-IN0EX2 

0020  GO   TC   401 

0021  402    FKEY<21=FSEY( 1)*FREY(2I 

0022  FCT(2)=FCT(1)+FCTI2» 

0023  INDEX1=INDEX1+1 

0024  NNN=N-INDEX1-IN0EX2 
00<;5  DL    -rUi     ill=l,ftrm 

0026  FCTIMIJ=FCT(IUH) 

0027  403    FREYIII I )  =  FREY 1 111*1 1 
0026  401    CONTINUE 

0  029  NN-N-INDEX1-INDEX2 

0030  CO   405    1*1 ,NN 

0031  CH12=(FCTI  II-FREYd  J)**2/FCT(II 

0032  CHI2=CHI2*J 

0033  405    CHI=CHI+CHI2 

0034  IFdNOEXl    .GE.    1  )IN0EX1  =  I  N0EX1*1 

0035  IFI1NDEX2    .GE.    1  U.\0EX2=INDEX2+1 

0036  DF=tJN-l-ESTPAR 

0037  RETURN 
003S  END 
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0001 

0002 
0003 
0004 
0005 
0006 
0007 
0008 
0009 
0010 
0011 


0012 
0013 
0014 
0015 
0016 
0017 
0018 

0019 
0020 
0021 

0022 
0023 

:::-. 

0025 
0026 
0027 
0028 

0029 
0030 
0031 


0032 
0033 
0  034 
0035 


SUBROUTINE    WBLFIT(MEAN,SI02,C.K,IER1 
C»»»    SUBROUTINE    COMPUTES    STARTING    POINTS    FOR    KUELLER'S    ITERATION 
REAL*4    K,MEAN 
CCMMCN/KBL/ALPHA 
EXTERNAL   F 
IER=0 

ALPHA*  1.0+SIG2/MEAN**2 

IF    (ALPHA. GT. 0.01    GO    TO    12 
PRINT    13,    ALPHA 
13    FORMAT! ■    ALPHA    IS    NEGATI  VE>'  ,G15.  7,  ■         l>C    SOLUTION    PCSSIELE'l 
IER»4 
RETURN 

C»*«    CALCULATE    INITIAL   STARTING    POINTS    WHICH    BRACKET    SOLUTION    FOR    K 
C«»«    SEARCH   FOR    THE    LEFT    HAND    STARTING    POINT 
12    00    10    1-1,100,2 

X=l.oool/I 

AA*F(XI 

IF(AA.GT.O.O)    GO    TO    15 

10  CONTINUE 

PRINT    11, X.AA, ALPHA 

11  FORMAT!'     CCULO    NOT    FINO    LEFT    HAND    STARTING    POINT    OR    RAN    OUT    OF    IT£ 
1RATI0NS',/'    X«',G12.4i'     F(XI=',G12.4,'       ALPhA= • ,G12.4I 

IER»3 
RETURN 
15    Xl-X 
C***    SEARCH   FOR    RIGHT    HANO    STARTING   POINT 
00    20    1*1,50 
X-XL+I 
>'.-F!X) 
IF    (AA.LT.O.O)    GO   TO   25 

20  CONTINUE 

PRINT    21, X.AA, ALPHA 

21  FORMAT!'    COUCO    NOT    FINO    RIGHT    HANO    STARTING    POINT    OR    RAN    OUT    OF    IT 
1ERATICNS',/'    X=',G12.4,'    F( X) ■' ,G12.4, '       ALPHA"  ,G12.4) 

I£R»3 

RETURN 
25  XR»X 
C 
C*»*  SOLVE  FOR  C  ANO  K 

CALL  RTMI(K,FN,F,XL,XR,O.COQOOI,100,IERI 

C=MEAN/GAMMA( 1.0+1.0/X) 

RETURN  ; 

ENO 
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0001  FUNCTION    WBCV) 

C*»*  SUBROUTINE  CALCULATES  VALUES  OF  EITHER  ( V/VRATEC)**»3*WSIBJLL  OR  JUST 
C     V.E1BULL  DISTRIBUTION 

0002  C0MMCN/LINK1/AIJ,BIJ,AA,BE,VMAXU»VRAT€C,FCTR 

0003  XN*tV/BU)**AU 

000*  IFIX8    -ST.    170.)    GO    TO   370 

0005  Ke=IAIJ/eiJ)»lV/BIJ)**<AIJ-ll*EXP<-XN) 

0006  RETURN 

0007  ENTRY    V3V.BCV1 

0008  XN=(V/BIJ)**AU 

0009  IFtXN  .GT.  170.)  GO  TO  370 

0010  V3WB»IV/VKATSD}*»3»UIJ/aiJ)*(V/BIJ)**( A  I  J- I ) «EXP(-XN) 

0011  GC  TO  371 

0012  370  *B=0.0 

0013  V3WB=0.0 

0014  371  RETURN 

0015  END 


0001  REAL    FUNCTION   F*4(XI 

0002  COMMON/ WBL/ ALPHA 

0003  REAL**    FF.XX 

C**»    EXTERNAL    FUNCTION    NEEDED    BY    SUBROUTINE    RTNl 

000*  xx=x 

0  005  FF=GAKMAt 1.0+2.0/XX) / I G AMMA I 1.C+ 1 . U/XXJ 1  "* 

0006  F«FF-ALPHA 

0007  RETURN 

0008  ENO 


FORTRAN  IV  G  LEVEL   21 


DATE  =  78135 


SUBROUTINE  RTH1 

PURPOSE 

TO    SOLVE    GENERAL    NONLINEAR    EQUATIONS    OF    THE    FORM   FCTIXI-0 
BY    MEANS    OF    MUELLER'S    ITERATION   METHOD. 

USAGE 

CALL    RTMI     (X,F,FCT,XLI,XRI,EPS,IENO,IERl 

PARAMETER    FCT    REQUIRES    AN    EXTERNAL    STATEMENT. 


10 
20 
30 
40 
50 
60 
70 
80 
90 


DESCRIPTION    OF    PARAMETERS 

X  -    RESULTANT    ROCT    CF    EQUATION    FCTIXI-O. 

F  -   RESULTANT    FUNCTION    VALUE    AT    ROOT   X. 

FCT  -    NAME    OF    THE    EXTERNAL    FUNCTION    SL3PRCGRAM    USED. 

XLI  -    INPUT    VALUE    WHICH    SPECIFIES    THE    INITIAL    LEFT 

OF    THE    ROOT 


RTH1 
RTMI 
RTMI 
RTHI 
RTMI 
RTMI 
RTMI 
RTMI 
RTHI 

RTHI  100 
RTMI  110 
RTMI  120 
RTMI  130 
RTMI  140 
RTMI  150 
RTMI  160 
RTHI  170 
BOUND  RTMI  180 
RTMI    190 


XRI 


I  END 
IER 


-  INPUT    VALUE    WHICH    SPECIFIES    THE    INITIAL    RIGHT    BOUNDRTMI  200 
OF    THE    ROOT    X.                                                                                              "MI  210 

-  INPUT    VALUE    WHICH    SPECIFIES    THE    UPPER    BOUND    CF    THE    RTHI  220 
ERROR    CF    RESULT    X.                                                                                    R™I  230 

-  MAXIMUM    NUMBER    OF    ITERATION    STEPS    SPECIFIED.  RTMI  240 

-  RESULTANT    ERROR    PARAMETER    CODED    AS    FCLLOWS  RTHI  250 

IER=0   -    NC    ERROR,  RTTII  260 

I£R*1    -    NC    CONVERGENCE    AFTER    I  END    ITERATION    STEPS    RTMI  270 

FOLLOWED    BY    IEHD    SUCCESSIVE    STEPS   OF  RTHI  280 

BISECTION, 
TFP  =  ?    -    BASIC    ASSUMPTION    FC T ( X L  1 1 »FCT I XR 1 1    LESS 

ThAN    CR    ECUAL    TO    ZERO    IS    NOT    SATISFIED. 


RTMI  290 

RTMI  300 

RTMI  310 

RTMI  320 

RTMI  330 

RTMI  340 

RTHI  350 

THERTMI  360 

RTHI  370 

RTMI  380 

RTMI  390 

RTMI  400 

RTMI  410 

RTMI  420 

METHOD                                                                                                                                                 RTMI  *3° 

SOLUTION    OF    EQUATION    FCT(XI=0    IS    CCNE    BY    HEANS    CF    MUELLER'S   RTMI  440 

ITERATION    METHOD    OF    SUCCESSIVE    BISECTIONS    AND    INVERSE                   RTHI  450 

PARABOLIC     INTERPOLATION,    WHICH    STARTS    AT    THE    INITIAL    BOUNDS    RTMI  460 

XLI    AND    XRI.    CONVERGENCE    IS    5UADRAT IC    IF    THE    DERIVATIVE    OF      RTMI  470 

FCTIXI    AT    ROOT    X    IS    NCT    EQUAL    TC    2ERC.    ONE    ITERATION    STEP         RTMI  480 

REQUIRES    TWO    EVALUATIONS    OF    FCTIXI.    FOR    TEST    ON    SATI SF ACTORYRTMI  490 

ACCURACY    SEE    FORMULAE    (3,41    OF    MATHEMATICAL    DESCRIPTION.           RTMI  500 

FOR    REFERENCE,     SEE    G.     K.    KRISTIANSEN,    2EPC    CF    ARBITRARY              RTMI  510 

BIT,    VOL.    3    119631,    PP. 205-206.                                                 RTMI  520 

RTHI  530 
'. RTMI  540 


REMARKS 

THE  PRCCECURE  ASSUMES  THAT  FUNCTION  VALUES  AT  INITIAL 
BOUNDS  XLI  AND  XRI  HAVE  NOT  THE  SAfE  SIGN.  IF  THIS  BASIC 
ASSUKPT1CN  IS  NCT  SATISFIED  SY  INPUT  VALUES  XLI  AND  XRI, 
PROCEDURE  IS  BYPASSED  AND  GIVES  THE  ERROR  MESSAGE  IER=2. 

SUBROUTINES  AND  FUNCTION  SUBPROGRAHS  REQUIRED 

THE  EXTERNAL  FUNCTION  SUBPROGRAM  FCTIXI  MUST  BE  FURNISKE 
BY  THE  USER. 


FUNCTION, 


SUBROUTINE  RTHI I  X ,F .FCT  ,XLI  ,XRI .EPS ,1  END, IER) 


RTHI  550 
RTMI  560 
RTMI  570 
RTMI  580 
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C 

PREPARE    ITERATION 

RTMI    590 

0002 

1ER="0 

RTHI    600 

0003 

XL-XLI 

RTMI    610 

0004 

XR=XRI 

RTNI    620 

0005 

X=XL 

RTMI    630 

0  006 

TOL  =  X 

RTMI    640 

0007 

F=FCTITOLI 

RTMI    650 

OOOS 

IF(F)lil6,l 

RTMI    660 

0001 

1 

FL-F 

RTMI    670 

0010 

X=-XR 

RTMI    680 

0011 

TCL=X 

RTMI    690 

0012 

F»FCT(TOLI 

RTMI    700 

0013 

■'■IF  12,16,2 

RTHI    710 

0014 

2 

FR»F 

RTMI    720 

0015 

c 

IF(SIGN(1.,FD-»SIGNI1.,FRII25,3,25 

RTMI    730 
RTMI    740 

c 

eASIC    ASSUMPTION    FL*FR    LESS    THAN    0    IS    SATISFIEC. 

RTMI    750 

c 

GENERATE    TOLERANCE    FOR    FUNCTION    VALUES. 

RTHI    760 

0016 

3 

1-0 

RTHI    770 

0017 

c 

c 

TOLF»100.»EPS 

RTHI    780 
RTHI    790 
RTHI    800 

c 

START    ITERATION    LOOP 

RTHI    810 

ooia 

c 

4 

1-1*1 

RTMI    820 
RTHI    830 

c 

START    BISECTION    LOOP 

RTMI    840 

0019 

DC    13   K«1,1EN0 

RTHI    850 

0020 

X=.5»CXL«XRI 

RTMI    860 

0021 

TOL»X 

RTMI    8f0 

0022 

F»FCT(TOLI 

RTHI    880 

0  023 

IF(FI5,16,5 

RTMI    890 

0024 

c 

5 

IFISIGNI1.,FI+SIGN(1.,FR))T,6,7 

RTHI    900 
RTHI    910 

c 

INTERCHANGE    XL    ANO    XR    IN    CROER    TO    GET    THE    SAKE    SIGN    IN 

F    AND    FR         RTMI    920 

0025 

6 

TOL=XL 

RTMI    930 

0026 

XL=>XR 

RTHI    940 

0027 

XR»TOl 

RTHI    950 

0028 

TOL'FL 

RTHI    960 

0029 

FL=FR 

RTMI    970 

0030 

FR»TCL 

RTHI    980 

0031 

7 

TOL»F-Fl 

RTHI    990 

0032 

A=F*TGL 

RTMI1000 

0033 

A=A*A 

RTHI1010 

0034 

IFtA-FRMFR-FL) 18,9,9 

RTMI1020 

0035 

S 

1F(I-IENDI17,17>9 

RTHI1030 

0036 

9 

XR«X 

RTMU040 

0037 

c 

FR-F 

RTHI1050 
RTMI1060 

c 

TEST    ON    SATISFACTORY    ACCURACY    IN    BISECTION    LOOP 

RTMU070 

0038 

TCL=EPS 

RTHI1080 

0039 

A.ABSIXRI 

RTHII090 

0040 

IFIA-1.  111,11,10 

R7NI1100 

0041 

10 

TCL=TOL»A 

RTHU110 

0042 

11 

IF(ABS(XR-XL>-TOL 112,12,13 

RTMIU20 

0043 

12 

IF(AES(FR-FL)-T0LF)14,14,13 

RTMU130 

0044 

13 

CONTINUE 

RTHU140 

c 

END    OF    BISECTION   LOOP 

RTHU150 

c 

RTM1U60 
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FORTRAN    IV    G   LEVEL 

21                                            RTH1                                     DATE    '   78135                        23/28/52 

C 

NO    CONVERGENCE    AFTER    IEN0    ITERATION    STEPS  FCLLCWED    EY    I  END                     RTMIH70 

c 

SUCCESSIVE    STEPS    OF    BISECTION   OR    STEADILY    INCREASING    FUNCTION             RTMI1180 

c 

VALUES    AT    RIGHT    BOUNDS.    ERRCR    RETURN.                                                                          RTMI1190 

0045 

IER.l                                                                                                                                                           RTNI1200 

0046 

14 

IF  (ABSIFRI-ABStFDI  16,16,15                                                                                                   RTMI1210 

0  047 

15 

X-XL                                                                                                                                                             RTHU220 

0046 

F.FL                                                                                                                                                             RTHI1230 

0049 

16 

RETURN                                                                                                                                                        RTMI1240 

c 

RTM1250 

c 

COMPUTATION    OF    ITERATED    X-VALUE    BY    INVERSE    PARASOLIC    INTERPOLATIONS TMI1260 

0050 

17 

A=FR-F                                                                                                                                                        RTNI1270 

0051 

DX«[X-XL]*FL»II.»F»IA-TCL)/[A»IFR-FLIII/TQ.                                                           RTNI1280 

0052 

XN-X                                                                                                                                                             RTMI1290 

0053 

FH»F                                                                                                                                                              RTHI1300 

0054 

X.XL-OX                                                                                                                                                      RTMU310 

0055 

TOL-X                                                                                                                                                           RTNU320 

005b 

F-FCTITOLI                                                                                                                                               RTNI1330 

0057 

IFIFI18,16,18                                                                                                                                       RTMI1340 

c 

RTHI1350 

c 

TEST    CN    SATISFACTORY    ACCURACY    IN    ITERATION    LOOP                                                 RTNI1360 

0058 

18 

TOL-EPS                                                                                                                                                      RTHU370 

0059 

A-ABS(X)                                                                                                                                                    RTMI1380 

0060 

IFIA-1. 120,20,19                                                                                                                               RTHI1390 

0061 

19 

T0L>70L»A                                                                                                                                                 RTM 11400 

0062 

20 

IF(AESICXI-T0L)21,21,22                                                                                                             RTHU410 

0063 

21 

IF(ABSIFI-T0LF)16,16,22                                                                                                             RTHI1420 

c 

RTKI1430 

c 

PREPARATION    OF    NEXT    EISECTICN   LOOP                                                                                 RTMU440 

0064 

22 

IF;ClCN!l.,r;:SI0N:i.,rL!I24,23,24                                                                                 RTMI1150 

0065 

23 

XR.X                                                                                                                                                              RTHI1460 

0066 

FR»F                                                                                                                                                              RTHI1470 

0067 

GO    TO    4                                                                                                                                                      RTHI1480 

0068 

24 

Xl-X                                                                                                                                                             RTHI1490 

0  069 

FL"F                                                                                                                                                             RTHU500                     i 

0070 

XR-XN                                                                                                                                                           RTNU510 

0071 

FR-FN                                                                                                                                                           RTBII520 

0072 

GO    TO   4                                                                                                                                                      RTMI1530 

c 

END    OF    ITERATION   LOOP                                                                                                                   RTHI1540 

c 

RTHI1550 

c 

RTHU560 

c 

ERROR    RETURN    IN    CASE    OF    WRONG    INPUT    DATA                                                                  RTK11570 

0073 

25 

IER=2                                                                                                                                                           RTHI1580 

0074 

RETURN                                                                                                                                                        RTBI1590 

O075 

ENO                                                                                                                                                                RTHI1600 
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APPENDIX  B 

2 
X   and  Power  Ratio  Tables 

2 
This  appendix  lists  the  remainder  of  the  tables  showing  the  % 

and  power  ratio  statistics  for  the  remaining  16  locations.   The  table 

format  is  identical  to  the  format  used  in  Table  2.6-2. 
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APPENDIX  C 
The  Program  BLOHARD 

The  purpose  of  this  program  is  to  find  the  economically  optimum 
WTGS  to  serve  a  particular  demand  load  given  sufficient  wind  speed 
data.   There  are  two  versions  of  BLOHARD  which  differ  only  in  the 
input  wind  speed  data.   Version  A  uses  the  observed  wind  speed  data  in 
its  calculations,  while  version  B  assumes  the  wind  speed  distributions 
are  given  by  a  beta  distribution.   The  program  is  written  in  FORTRAN  IV 
for  use  on  the  Kansas  State  University  ITEL  AS/5  System  (equivalent 
operationally  to  an  IBM  370/158) .   Because  of  the  liberal  amount  of 
comment  cards  in  the  listing  which  follows  and  the  variable  names  of 
high  mnemoic  content,  only  a  brief  modular  description  of  the  program 
is  given.   Following  the  program  listing  is  a  sample  output. 

BLOHARD  maximizes  the  objective  function  of  Eq.  (3.2-22)  according 
to  the  methodology  described  in  Section  3.2.   The  subroutine  OBJN  per- 
forms the  calculations  necessary  for  the  optimization  methodology.   The 
WTGS  cost  model  is  calculated  in  the  subroutine  MONEY.   Although  only 
one  cost  function  is  specified,  MONEY  is  capable  of  using  another  cost 
function  of  the  same  form  as  the  functions  given  by  Eqs.  (3.2-20)  and 
(3.2-21).   The  calculations  required  in  the  simplex  search  pattern  are 
performed  by  the  subroutine  SIMPLEX.   Finally,  the  BREZE  subroutine 
determines  the  maximum  and  minimum  power  requirements  and  maximum 
wind  speed  of  the  input  data. 
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FORTRAN  IV  G  LEVEL   21  HAIN               DATE  -  78135         23/29/01 

C*******»t*w*»***«t«*»t****>*i»i*  Bl CHARD  *•»«»**»***'*•*****»*»*«***•**»#*»♦* 

c» 

C»  THIS    PROGS AN    MATCHES    DAILY    HINT    SPEED    01 STRi BUTIChS    TC    THE    REOUIREO    DAILY 

C»  LOAO    DEMAND.       THE    WTGS    IS    OPTIMIZED    TG    FIND    THE    BEST    SUED    SYSTEMU.E. 

C*  OPTIMUM    RATED    PD.EP.    AND    RATEC    SPEED)    SO    AS    TC    MAXIHI2E    THE    ELECTRICAL 

C»  SAVINGS.      BETA    DISTRIBUTION    OR    OBSERVED    HIND    SPEED    DISTRIBUTION    CAN    BE   USEO 

C» 

C« 

C»  INPUT    DATA:        (VERSION    A   -    EETA    DISTRIBUTION    AS    KIND    SPEED    MCDEL) 

C« 

C»  CAPO    I           FORHAT    (12) 

C»  NHALF    .    THE    HALF    VALUE    OF    THE    EVEN   ORDER    GAL'SS-LEGENDRE    OUADRATURE 

C*  USED    TO    EVALUATE    THE    NECESSARY    INTEGRALS. 

C» 

C*  CARC    2           FORHAT    (4920*01 

C*  ROOTMi     -    CUADRATURE    ORDINATESIONLY    POSITIVE    VALUESI 

C»  (HAY    BE    MANY    CAROS) 

C» 

C*  CARD    3           FORHAT    (AG20.01 

C»  HEIGHT!  II    =    CUADRATURE    HEIGHTS 

C*  (HAY    BE    MANY    CAROS I 

C» 

C*  CARD  A            FORMAT     (£I5,2G10.0I 

C»  NDINT    •    NUMBER    CF    DAILY    SUBINTERVALS    FCR    HUE    AND    LOAD   DATA 

C*  NYINT    a    NUMBER    CF    SEASONS    FOR    WHICH    HIND    AND    LEAD    DATA   ARE    GIVEN 

C*  NYEARS    *    NUMBER    CF    YEARS    HTGS    IS    AMORTIZES 

C*  INT    =    YEARLY    INTEREST    RATE 

C»  ICCST    «    NUMBER    CF    WIGS    COST    FUNCTION    LSEOCTkO   ARE    POSSIBLE! 

C*  VREFU1     t   REFERENCE    RATED    SPEED    HF    FIRST    WTfS    CrKT    FUNCTIONIKNOTSI 

C*  VREFI2)    ■   REFERENCE    RATEC    SPEED    CF    SECOND    HTGS   CCST    FUNCTICNIKNOTS J 

C» 

C*  CARD    5           FORMAT    IBG10.01 

C»  2(1)    -    INITIAL    SIMPLEX    POINT    RATED    PCWER(KU) 

C*  Z(2)    =    INITIAL    SIMPLEX    PCINT    RATED    SPEED(KNCTS) 

C»  EPSI    «    CONVERGENCE    CRITERION    OP    SIMPLEX    TECFNICUE 

L*  FRACT1    =    FRACTION    OP    LOAD    MAXIMUM    PDkER    DEMAND    USEO    AS    STEP    SUE    IN 

C»  SIMPLEX    TEChNIQUE 

C»  FRACT2    =    FRACTION    OF    MAXIMUM    HIND    SPEED    USED    AS    STEP    SIZE    IN    SIMPLEX 

C»  TECHNIQUE 

C* 

C»  CARD    6           FORMAT    (20AAI 

C«  TITLE    -    TITLE    CARD    FCR    PROBLEM   TO    BE    ANALYZED 

C» 

C*  CARD    7           FORMAT     (8G10.0) 

C»  MEAN    »    MEAN    HIND    SPEEO    IN    M-TH    DAILY    INTERVAL    OF   MM-TH    SE4S0NIKN0TSI 

C*  VAR   «    VARIANCE    IN    HIND    SPEEDS    IN    M-TH    DAILY    INTERVAL    CF    MM-Th    SEASON 

C»  [KU!ITS**21 

C»  VMAX(K,MM1    •    HAXIMJM    HIND    SPEED    IN    M-TH    DAILY    INTERVAL    OF    HM-TH    SEASON' 

C*  (KNGTSI 

C*  (MAY    BE    MANY    CAROS) 

C» 

C*  CARD   8           FORMAT    (8010.41 

C»  PLCAOIl.JI    =    AVERAGE    POWER    REQUIREMENT    FOR    THE    I-TH    DAILY    INTERVAL    IN 

C«  J-TH    SEASCN(tCk).       (HAY    EE    MANY    CARDS) 

C* 

C*  CARD    9           FORMAT    (2G10.0) 

C*  SBUY   «    COST    OF    PURCHASED    ELECTRICI TY I i/KWHI 

C»  tSEL    '    VALUE    CF    ELECTRICITY    SOLDIJ/KWHI 
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C» 

C» 

C» 

INPUT  DATA 

C» 

C» 

CARD  1 

c» 

INN 

c» 

c» 

CARD  2 

c> 

ivin 

c« 

c» 

CARO  3 

c» 

c» 

c* 

CARD  A 

c« 

c« 

c» 

CARD  5 

c» 

c* 

HAIN 


DATE 


78135 


tVERSION  6  -  OBSERVED  HINO  SPEED  CISTRIEUTION  USED) 

SUBINTERVALSIUSUALLY  11') 


FORMAT    1131 

TOTAL    NUMBER   OF    POSSIBLE    SPEE 


FORMAT    (12F5.2) 
(II    =    ENDPOINTS    OF    WIND    SPEED    SUB  INTERVALSIKNCTS) 


FORMAT    I5I5.ZG10.SI 

(SAME    AS    CARO    A  -    VERSION   Al 


FORMAT  (8G10.0I 

(SAME  AS  CARD  5  -  VERSION  Al 


FORMAT  (20A4) 
[SAME  AS  CARD  6 


VERSION    A) 


CARD    6  FORMAT    1 12, 13, 1014  , 1 E, IS! 

MCNTH   =    MONTH    FROM    WHICH    SPEED    DATA    IS   OBTAINED 

NTIME    =    DAILY    TiME    PERIOD    FROM    WHICH    WIND    SPEED    CATA    IS    OBTAINED 

IFREQ(I)    =    FRECUENCYU    1C00I     OF    CBSERVATICNS    IN    l-TH    SPEED    SU5INTERVAL* 

(BEGINNING    WITH    FREQUENCY    IN    2ND    SPEED    SUB1NTERVAL1 
NUKMY    =    SPACE    FOR    DATA    IDENTIFICATION    PURPCSES1CAN   ALSO    BE   USEC   TO   ACO» 

TWO   MORE    SPEED    SU3 INTERVALS) 
NCES    ■    TOTAL    NUMEER    OF    WING   SPEED    DATA    OBSERVATIONS 
(MAY    BE    MANY   CARDS) 


0001 
0002 

0003 
000* 
0005 
0006 

0007 

oooe 
ooos 

0010 


c* 

(SAME  AS  CARD  3 

•  VERSION  A) 

C» 

c« 

CARO  8 

FORMAT  (2G10.0) 

C*  ISAME    AS    CARD    9    -    VERSION    A) 

C* 

C*    WRITTEN    BY    L.    A.    PCCH ,    KANSAS    STATE    UNIVERSITY,    JANUARY    1978 

C* 

£*=***  ***********************************************  ft*************************: 

C 

IMPLICIT    REAL»8(A-H,0-Z,$) 

P.SAL*  S    FRECK1)  ,  VI 4 1  ) ,  FIS1 )  ,T  ITLEI20I  ,  I  VINT  (41 J , 
X    FCTM6»4)  .A(8,4)tS(8»4),Z<2>  .STEP  (  2 )  ,  PLCAO  (  5  ,  <■  J  ,V:J.AX(8,4) 

REAL*8    Y(41),RC0T(20)  ,WEI  GHlf  20  1  ,X(  41  J  ,HN(6  ,i)  ,  VR (  3 , 4) ,  VREF  (2) 

REAL*6    HKEAf,,HAXSAV|MEAN»K 

INTEGER**    IFREQI21) 

CQMMCN/i.INKl/A,S.PL0AD|FCTR(VMAX,VSEF,SDYS|HRS,*3L,Yf  SSELtCRFf 
1    VCUTIN,NYINT,NDINT,ICOST 

CCKMCN/LINK2/AIJ,£lJ,VMAXIJ,VRATED,PfiATED,FRACT3(FRACT4 

CCMMCN/LINK3/RCCT,WEIGHT,NHALF 

REAL* 5    FRECNf  12,8,4)  ,  VI  NT ( 12 , 8 ,41  ,MIM<8,4) 

C0MMCN/LINK4/FREQN,V1NT,MIN 
C***    REMOVE    COMMENT    IF    VERSION    A    IS    TO    fiE   USEC    ANC    ACC    COMMENT    TO    APPROPRIATE 
C  VERSION    B    CARDS 

C  READ15,113JNHALF 

C     113    F0RMATU21 

C  READ  (5,  12HR00TII),I  =  1,NHALFI 

C  READ  I  5, 12  J (HEIGHT! I )  ,  I=l,NHALFI 

C       12    FQRHATI4G20.0) 
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0011 
0012 
0013 
0014 
0015 
0016 
0017 
0016 
001? 
0020 
0021 
0022 
0023 
0024 
0025 
0026 
0027 
0028 
0029 
0030 
0031 
0032 
0033 
0  034 
0035 
0036 
0037 
0C3S 
0039 
0040 
0041 
0042 
0043 
0044 
0  045 
0046 
0047 
004» 


C 

c*** 

c 

c 

c 

c 


C  READ    10,NDlNT,NYH,T,NYEARS,INT,ICCST,VREF(ll,VREFt2J 

C       10    FCRMAT15I5.2G10.0I 

C  READ    l,Z(lliZI2i ,ERSI.FRACT1,FRACT2 

C    111    READI5.100I    TITLE 

C    100    FCRMAT(20A4> 

C  kRITE{6,110)    TITLE 

C    110    FCRMATl  '1' ,20A4I 

C  DC    170    MM=1,NYINT 

C  DO    170   H-l.NDINT 

C  READ1,MEAN,VAR,VHAX(H,MN) 

1  FORMATI8G10.0I 
FIT    DATA    TO    A    BETA   DISTRIBUTION    BY    HATCHING    MEAN    AND    VARIANCE 

A[M,HM)=(MEAN/VMAX(M,MM1)*IMEAN*(VMAX(M,MM)-MEAN)/VAR-1.0DO) 

6(M,MK)=f  VMAXIM,MV1-MEAN)«AIM,MM)/MEAN 

VVAR=AtM,MK>»BlK,M-l)»VKAX(M,MMI»«2/[(A(«,MKl»Bll'.,HM>l»»2»<A<H,MM) 
X    *B(F,HH)+I. 000)1 

MMEAN=A(M,MM)*VMAX(H,MM)/  [AM,MM]+8IM,HK)) 

FCTR(M,MM]  =  DGAM,".A(A(M,MH)-»B(M,MH)  1  /  (  VMAX(  M,  MH1  »DGAMMA(  A I  H.HNI  ) 
X  »DG-AKMA(BtH,HH) )) 
***  BEGIN  VERSION  B 

READ(5,2)  INN 

2  FCRMATII3I 
NN-INN+1 

READ!  5,111  [IVINTU),I-1,NM 
11    F0RHATI12F5.21 

READ    10. KD1 NT, NY  INT, NYE ARE, INT, ICDST.VREFUJ ,VREF{2) 
10    FCRKAT(5I5,2G10.0I 

READ    1,2(  11,21 21, EPSI,FRACT1,FRACT2 
1    FCRMAT'.8GI0;01 
111    READ(5,1001    TITLE 

100  FCRMAT(20A4) 
kRITE(t,110>    TITLE 

110    FCRMATl '1' ,20A41 

OC    170    MK=1,NY1NT 

DO    170   M-l.NOINT 
99    N-INN 

DC    172    I-l.NN 
172    VINT<I,M,MKI«IVINTUI 

VMAXtF.MHI* IVINT(NN) 

READ!5,101>MGNTH,N7IM£,(IFRE0(I),I-2,NI,NUMMY,    NOBS 

101  FCRMATl 12,12,1014,18,151 
DC  401    1-2, N 

401    FREQIII-OFLOATUFREQIIH/IOOO.OO 

SUM-O.OOO 

DC   54    1-2, N 
54    SL'M-SUMIFREOII) 

FP£0(1)-1.0D0-SUH 

SUH-SUH»FREC(1> 

SUMI-O.ODO 

DC   41    1-1, N 

5UMl=SUMl»FRE0tI> 
41    Ftll-SOMI 

OC    130    I-l.N 

IF(FII)    .GE. 

130  CCNTINUE 

131  N-I 
MIM<M,MM)*N 
SUM2-0.000 


9.9999990-01)    GO    TO    131 
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0049 
O050 
0051 
0052 
0053 
0054 


0055 
0C5t 
0057 
0058 
0059 
0040 
0061 
0062 

0063 


0  064 
0065 
0066 
0067 
0066 
0069 
0070 
0071 

0  072 
0073 
0074 
0075 

0076 
0077 
0078 


0079 
0080 
0081 
0082 
0083 
0064 
0085 
0086 
00S7 
0C88 

0089 
0090 
0091 
0092 
0093 
0  094 
0095 


CO   46    IK=1,N 

SUM2=SUM2+FREQtIKl 

FI1K]»SUN2 

Y!  1K)=V!NT(  1K+1  ,M,MK1-VINT(  IK,H,MN> 

V<m=C.5D0*(VINT(IX*l,M,PMl+VINTUK,N,KK)) 
46    FREQN[IK,M,HMI  =  FREQII!0/YC1K» 
C 
C*»*    CALCULATE    MEAN    AND    VARIANCE    OF    HIND    DATA 

MEAN=0.0DO 

VAR=0.000 

CO    30    I'liN 

30  NEAN=HEAN+V11>*FRE0I II 
DO    31    I-l.N 

31  VAR=VAR    +    (V(I>-M£AN)**2*FREQIII 
MNIM,fMI=H£AN 

VPtM^HI-VAR 
C***    END   VERSION    B 

170    CONTINUE 
C 
C*»*    READ    IN    LOAO    DATA 

HRS=24.0DO/NDINT 

SDYS=365.OD0/NYINT 

SHRS=24.0D0*S0YS 

00    21    J*1,NYINT 

21  READ    112,     CPLOACIl,J],I=liNDINTJ 
112    FCRMAT18G10.4I 

PRINT    22.HRS 

22  F0RMATI//'0AVERAGE    VELOCITY    FOR    ENTERPRISE    DURING    EACH',G10.3, 
1     'HOUR    INTERVAL    BEGINNING    AT    MIDNIGHT    (IN    KN0TS1:'! 

00    26   J=1,NYINT 
26    PRINT    16,J,IKN(I,J>  ,I»l.NUNTI 
PRINT    23,    HRS 

23  FORMAT!  // ' CELECTRICAL    POKER    DEMAND    FCfi    ENTERPRISE    DURING    EACH' 
1,G10.3,'HOUR    INTERVAL    BEGINNING    AT    MICNI.hT    UN    KWI : "  I 

DO    25   J=1,NYINT 
25    PRINT    16,J,(FLCAD(1,J>,I=1,NDINT> 

It   FORMAT! ■    SEASON', 12,':         • , 10G1 1.4 ,/( 13X.10G11 .4) I 
C 

C**»    CALCULATE    AVERAGE    POWER    NEEDS 
HAV=0.0DO 
PAV=O.ODO 
DC    300    J-l.NYINT 
DO    300    I-l.NDINT 
WAV=WAV+MN( I,  J) 
300    PAV=PAV+    PL0A0U,JI 
WAV=WAV/(NDINT*NYINTI 
PAV=PAV/(N01NT«NYINTI 
PRINT    311,P»V,kAV 
311    FORMAT  I //'OAVERAGE    POWER    REQUIREMENT* ■ ,G12. 5 , ■ KW.  ' ,5X, 
!     'AVERAGE    WINE    SPEED= ' , G12.5 , 'KNOTS ' I 
CALL    en  =  ZE(PLOAD,NaINT,NYINT,PMAX.,PMIN,VMAX,VMAXX) 
ppiNT    *4,     PMAX,PH1N 

34  FOHUTC    EXTREMA    CF    POWER    RECUI  REMENTS^  ARE' ,2G15.8,  "KW.  ■  I 
1111    REAOI5,3S0,EN0-98)     5BUY,»SEL 

350    FCRMATI2G10.0I 

PRINT    35,iBUY,$SEl 

35  FORMAT!//'     ASSUMED    COST    OF    ELECTRICITY   -     '.G10.3,'    [J/XWH1V 
1    '    VALUE    OF    ELECTRICITY    SCLD    -    '.G10.3,'     (S/KWM'I 
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C 

C**   CALCULATE    OPTIMUM    SAVINGS   AND    fcTGS    PARAMETERS 

0096  PRINT    33 

0097  33    FORMAT(//'OGEN.    RATING  VRATEO  VCUTIN  ANN.    WIND 

XCOST  GEN.    ELECT.  PURCH.    ELECT.  SCLO/UASTED  ANN.    NET    SA 

XVINGS'  ,/,5X,MKWI«  ,T18,  ■  (  KT  S I  '  ,T2E,  '(KTS)  ',  T53, »(il ■» 

X    T6S,'tKUH>  I  ,T65,'(KUH)"  ,7101,'  UWril'  ,T120."(1I'I 

0098  SPCAF=(1.DO»INT/1CO.DOI»»NYEARS 

0099  CRF=SPCAF»INT/100.DO/(SPCAF-1.DO> 

oioo  ziii-zin 

0101  2(21=212) 

0102  S7EPI  n=PMAX*FRACTl 

0103  S7EP [ 2)=VMAXX*FRACT2 

010*  CALL  SIMPXI Z,MAXSAV,2,STEP,1CC,100,  EPSI  ,1  .000,0. 5D0, 2. 0001 
C 

C».«  PRINT  FINAL  ANSWERS 

0105  VCUTIN=A. 6*15890-01*2(21 

0106  MAXSAV=-MAXSAV 

0107  PRINT    4*, 2(11, ZI2I, VCUTIN, MAXSAV 

0108  **    FCRMATCO' .//'ODPTIMUM    GENERATOR    SI  2E='  ,  G10  .3,  'KW* ,  SX,  'RATEO    SPEEO 

1  =  '.G10.3,UKTS)',5X,'CUT-IN    SPEEDs'.GlO.S.'IKTSI'.SXi 

2  /'  MAXIMUM  SAVINGS  ($}=*, G15. 8) 

0109  GO  TO  1111 

0110  98  S1CP 

0111  END 
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0002 
0003 

000* 

0005 

0006 
0007 
O008 
0009 
0010 
0011 
0012 
0013 
0014 
0015 
001b 
0017 


0016 
0019 
0020 
0021 
0022 
0023 
0024 
0025 
0026 


SUBROUTINE    OB JN[ X.NETSAV.NI 
>•    SUBROUTINE    COMPUTES    VALUES    NEEDED    IN    ECONOMIC    OBJECTIVE    FUNCTION 

IMPLICIT    REAL"8IA-H,C-Z,*I 

REAL*£    X(2),A[a,4),Bi6,41  »PIGA01S,*) ,  FCTP.  (8,4)  ,DEFICTIB,4), 
1    BUY(e,4),SELIE,4l,£XC£SS(8,4],GENCB,4l,VMAXl6,4I,VREF(2l 

R£AL»6    NETSAV, LIMIT 

CCMMCN/LINKl/A,B,PLOAD,FCTR,VMAX,VKEF,S&rS,HRS,»BUY,»SEL>CRFt 
1    VCUTIN,NYlNT,t.DINT,ICOST 

COMMOK/LINK2/AIJ,B1J,VMAX1J,VRAT£0,PRATED,FRACT3,FR«CT* 

REAL»6    FRECNI 12,8,41 ,VINT ( 12, 8, 41 ,MIM( 8,4) 

C0MM0N/LINK4/FRECN,VINT,MIM 

PRATED=X( 1) 

VRATED-XI2I 

IFIVRATEO    .GT.    O.DO    .AND.     PRATED    .GT.    O.DOI    GO   TO    10 

NETSAV  1.0010 

RETURN 
10    VCUTIN=4.641589D-01*VRATED 

SUMBUY=C.OD0 

SUMSEL=0.OD0 

SUMGEN=C.OD0 
«•    REMOVE    COMMENT    IF    VERSION    A    IS    TO    BE    USED    AND   ACC    COMHENT    TO    APPROPRIATE 
VERSION    B    CARDS 

EXTERNAL    FI.V3FI 

DC   40    I=1,NDINT 

CO   40   J*1,NYINT 

AIJ'AII.JI 

r? ij=bi I . j> 

VMAXIJ-VHAXII.JI 

1FIPRATED    .GT.    PlOADIl.J))    GO    TO    20 

BUY! I, J!         'PLCADII. J)-PRATEO»FCTRI I , Jl »1 GLCUACIVCUTIN, 0MIN1 
1    (VRATED,VMAXIJI,V3FIlfGLCUAD(VRATED,VMAXlJ,FII) 

SELIl.JI'O.ODO 

GO    TO    30 
20    VO=VRATED»(PLOAD 1 1,  J  I /PRATED) «»3. 3333333223333330-01 

LIM1T=DM!N1(V0,VMAXIJ) 

IFIVO    .LT.    VCUTIMLIMIT=DPAX1(V0,VCUTINI 

IFIVCUTIN    .GE-    VMAXIJI    LIMIT  =  VMAXU 

BUYII.JI         .FCTRU,JI»(PLCAC(I,J)»GLCL'«C!C.000,LIHIT,F1I 
1    -PRATED*GLCUAD(  VCUTIN.DMINU  VO  .VMAXU)  ,  V2F  1 1  ) 

SELII.JI         .FCTR(l,JI*lPRATED»IGLuU.\DlOMAXl(VC,VCUTINI,DMINl 

1  (VRATED.VMAX1J) ,V3FI) 

2  ,GL0UA3IVRATED,VBAXIJ,Fll>-PL0Aull  ,  JI«GLCUAD(  CM.AX1 1 VD,  VCUTINI  , 

3  VMAXIJ,FII) 

30    GENII, J)=PLCAD(I,J>-euYU,JI 

5UMBUY=SUMBUY+8UY(I ,J) 

SUKS£L"SUKSEL+SEL(I,J1 
40  SUMGEN=SUMGEN+GENU  ,J1 
**    BEGIN    VERSION    8 

DC    306    H=l ,NDINT 

CO    306    MM=1,NY!NT      ' 

MI=MIMH,MH) 

H3UY=O.0DO 

HSEL»0.000 

IFtPRATED    .GT.    PLCADIM.MM))     GG   TO  320 

CO    300    1-1 .HI 

IF(V  INTl 1  +  1,M,MM)     .LE.    VCUTIN1G0   TO    300 

1FIVINTC1 ,K,MM)     .GE.    VRATEOI     GO    TO    304 
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0029 

OO30 

0031 
0032 
0033 
0034 

0035 
0036 

0037 
0036 

0039 
0040 
0041 
0042 

0043 
0044 
0045 
0  046 
0047 
0048 
0049 
0050 
0051 

0052 

0053 

0054 

0055 

0056 

0057 

0058 

0059 
0060 
0061 
0062 


IV  G  LEVEL   21 


OBJN  DATE  -  7E135 

.LE.  VCUTIN  .ANO.  VINT!  J*li  N,MM  .GE. 

.LT.    VCUTIN    .AND.    V  IN7I 1*1, K,MK]    .LT. 

.LE. 

.07. 


23/29/01 

VRATEO) 
VRATEOI 
VRA7E0I 
VRATEOI 


IF(VIN7U,M,MMI 
1    GO    TO    302 

IFIVINTII  ,M,MM] 
1   GO    TO   303 

IFIVINTII, K.MHI     .GE.    VCUTIN    .AND.    VINTI  1*1. KtHp.l 
1    GC    TC    305 

IFIVINTII  ,M,MM)    .GT.    VCUTIN    .AND.     V1NTI  1*1  ,K,M* ) 
1G0   TG  333 

PRINT    32, i  ,M, MM, PRATED, VRATED 
32    FORMAT! 'OINTERVAL    ODES    NOT    FIT    ANY   CATEGORY' ,/ ,315,2G15.7J 

GC    TC   300 
333    HBUY=HBUY+( 0.2  5DO*(VR A7ED**4- VINT ( I,M, MM  1**41 /VRA7EC**3* 
1    VINT  (I+1,M,MM)-VRATED)*FREQNII  ,M,HHI 

GC    TC    300 
30  2    HBUY  =  H5UYM0.25D0*I VRATED**4-VCUTIN**4) /VRATED**3*VINT(  I+1,M,MN)- 
1    VRATED)*FREQN(I,M,HM 

GO    TO    300 

303  HBUY  =  HBUY-*.25DO*FREOHU,M,MK)*(VINT1I*1,N1MM)*»4-VCUTIN**4I/VRATEO 
1    **3 

GO   TO    300 

304  HBUY=KBUY*FRE0NI!,H,MM)*(VINT!H-1,M,MM)-VIN7II ,M,M1I 
GC    70   300 

335    HB>JY=HEUY  +  0.25D0*FREQN(I  ,  M  ,  MM)*  I  VINTI  1*1 ,*,  MM  1**4- 

1    VINTtl.MfMM)**4)/VRATED**3 
300    CONTINUE 

BUYIM,MH)=PLCADtM,MM)-HBUY«PRAT£D 

SELIK,MM)*HSEL 

GO    TO    330 

320  VD=VKAI  tU'l^LUA^lM.MMj/CRArtU)**^.  iiiiiiU-(il 
DO    301    1  =  1  ,MI 

IF(VINT[I+1,M,MM)    .LE.    VCUTIN)    GO    TO    311 

IFIVINTI  I,K,MM)     .GE.    VRATED)    GO    TO    318 

IFIVINTI I. M, MM)     .LE.    VCUTIN    .AND.    V  INTI I*1,H,MM)     .GE.   VO    .AND. 
1    VINT(!+1,M,MM)     .LE.    VRATED)    GC    TO    312 

IFIVINTII, M, MM)     .LT.    VCUTIN    .AND.    VJKTC 1*1, M, MM)    .LE.    VO) 
1    GO    TO    313 

IF(VINT(I,H,MM)     .GT.    VCUTIN    .AND.    VINTI  1+  l.M.MK)     .L£.    VO) 
1    GO    TC    314 

IFIV1NTI  I,M,KM)     .LE.    VD    .ANC.    VINTI  1*1, H,KM)    .LE.    VRATEDI 
1    GC    TC    319 

IFIVINTII ,M, MM)     .LT.    VCUTIN    .AND.    V  INTI  Hl.M.MF  I    .GT.    VO    .ANO. 
1    VINTU  +  1,M,MM]     .GT.    VRATED)    GO    TC    321 

IFIVINTII ,N, KM     .GT.    VO    .ANC,    V  INT  1 1*1, M, MM  I    .LE.    VRATEOI 
1    GO    TO    315 

IFIVINTII, H, MM)     .GT.    VD    .ANC.    V 1NTI 1*1, M, MM I    .GT.    VRATEDI 
1    GO    TO    316 

IFIVINTII,M,MH]     .L7.    VO    .AND.    VINTI 1*1 ,N, MM1     .GT.    VRATEOI 
1    GO   7C    317 

PRINT    31,1 ,M, MM, PRAT EO, VRATED, VD 
31    FCRMATI 'OINTERVAL    DOES    NOT    FIT    ANY   CATEGORY' ,/  ,31 5, 3G15.7I 

GC   TC   301 

321  HBUY=HBUY*FREONtl,M,MM)*(FLCADlM,MMJ*{DMA*l I VO, VCUTIN )- 

1  VI  NT  II ,M,MM))-PRAT£0*0.2SDC°IDMAX1[VD,VCUTIN)*»4-VCUTIN**4)/ 

2  VRATED**3) 

HS£L«=HSEL+FREON<  I  ,M  ,MM)  *(  PRATED*!  0. 25D0* I VRATED**4- 
lDMAXl!VD,VCUTINl*»4)/VRATED»*3+ViNTll*l  ,M ,MM1-VRATEC)-PL0A0(M,MMI* 
2    I VINTI I*1,M,MM)-DMAX1(  VD, VCUTIN) )) 

GC   TO   301 
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OBJN 
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23/29/01 


0065 
0  066 
0067 


0069 
O07O 


0071 
0072 


0073 
0  074 


0075 
0076 


0077 
0078 


0060 
0081 

0082 
0083 

0084 

OC85 
0086 
0087 
0086 
0089 
0090 
0091 

0092 
0093 
0094 
0095 
0096 
0097 
0098 
0099 
0100 
0101 
0102 


311  HBUY«HBUY»FREON(I,H,KM>*PLOAD(M,MM)*(VINT[I*1,I<,MMI-V1NT(1,M,HM1I 
GC    TO    301 

312  HBUY=HCUY+FREQNI  I,M,MM)»tFLCAClM,MM)*(OMAXl(VC,VCuflN>- 

1  VINT  II  ,M,HM) )-PRATED*C.2500*(DMAXl[VCUTIN,VD)**4-VCUTIh+*41/ 

2  VRATED**3) 

HSEL=HEEL  +  FREeNU,M,MM)*IPRATED*O.25D0*(  V1NTU-,1,K,HHI*«4- 

1  DMAX1(V0,VCUTINI**4)/VRAT£C»*3-PL0ACIM,MK]*( VINT (1*1, H,MHI- 

2  DMAX1I VC, VCUTIN))) 
GG  TG   301 

313  H6UY=H5UY+FRE0NII  ,M,MM)*{ FLGADI M, MM )* [V INT(  1+1 ,M, MM )-VINT( I ,M,MH> ) 
1    -PRATEC»0.25DC*(VINTU»l,M,MMJ**.,-VCUTIt>**4l/VRATED*»31 

GC   TC  301 

314  H8UY=HSUY+FRECN(I  ,M  f  MMI  *  (  FLCADt  H  ,  MM  )  *  1 V  INT  (  1+  1  ,  M  ,  MM)-V  I  NT  ( I  ,  M,  MN>  I 
1    -PRATEC*0.25D0*tVINTU-*l,M,MM)*»4-VINTII  ,H  ,MH)»*4I /VRATED**3) 

GO    TC   301 

315  HSEL=KSEL+FREQNII,M,MM)*(PRATED*0.25D0* I  VINT ( 1  +  1 ,H,MM)**4- 

1  vlNT(I,H,MMI*»4)/VRATEC*»3-PL0AD(M.,MM)*IVINTIl»l,NtMNl- 

2  V1NTU  ,M,MMI  I  I 
GC   TO    301 

316  HSEL=HSEL  +  FRECN(  I  ,M,M.MJ*t  PRATED*!  .25D0*(  VRATED**4-\'INT(  I,  M,  MM  1**41 

1  /VRATED** 3  +  VI NT C 1  + 1, M, MM  I- VRATED I-PLCAC I *, MM)* (VINT II*1,M,  HK)- 

2  VINTII,M,MM)I) 
GO   TC   301 

217    HBUY«H£UY*FREON(I ,M ,MM) »( PLCAD( M,"M)*(VD-VINT ( I ,M,MM) l-PRATEO* 
1    C. 2 500* IV 0**4- VI NT  I  I , M,HM) **4) /VR ATE 0**3  I 

HSEL=HSEL+FRECN( I . M, MM) *( PRATED* 1 .2  5  00* ( VRATE0**4-V 0**4 )/VRATE0«*3 
1    +VINTE  1+1,M,MM>-VRATED)-PLGAD(M,KH)*(V1NT[I+1  ,.".,  MM1-VD)  I 
GC   TC   30i 
213    HSEL=HSEL+FREQNl I ,M , MM) *l PRATED-PLCAOIH.HM )*IVINT(I+1,M,MMI- 
1    VINT(I,M,MM1I 
GO   TO   301 
319    H5UY=HEUY*FRE0N( I ,M , MM) *( PLGADI M, "Ml* t VC-VINT ( I ,N,MM) I -PRATED* 
1    Q.2SB&*(V0**4-VINT ( I,M,MM)*»4»/VRATEO*»3) 

KSEL-HSEL+FRECM1  ,  M,  KM)  »(  PRAT  ED*  0.2  5  DO*  I  VINT  I  1*1,  H,MM  )**4-VD**4)/ 
1    VRATED**3-PL0AD[M,MM)*(VINTtI*l,M,MM)-V0>) 
301    CONTINUE 

SEL(M,MM)=HSEL 
BUY{K,MM»*HBUY 
330    GENIM,MM]=PLCADIM,MM)-BUY(H,MM 
SUMEUY=3UMBLY+BUYIM,MHI 
SUMSEL=SUMSEL+SEL(M,MH1 
306    SUM&EN- SUMGEN*GEN(H,MMI 
***    END    VERSION    8 

SEAHRS=SDYS*HRS 
SUMBUY=SUMBUY*SEAHRS 
SUMSEL=SUMSEL*SEAHRS 
SL'MGEN  =  SUMGEN*SEAHRS 

CALL    MG;,EY(PRATED,VRATEC, VREF.      ANCC5T,CM,CRF, ICOST) 
NETSAV=-SUMGEN*SBUY  +  ANCCST-S'JMSEL*SSEL*CN 
ANS=-NETSAV  '■ 

PRINT   4 3. PR AT  ED, VR AT  ED, VCUTIN.ANCCST, SUMGEN.SUKBUY ,SUMSEl,ANS 
43    F0RHAT(2X,3G11.5,11X,G13.£, 5X,G13.6,4X,G13.6,2X,G13.6,5X,G13.61 
RETURN 
END 
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0002 
0003 
0004 
0005 
0006 
0007 
0008 
0009 
0010 
0011 
0012 
0013 
OOU 
0015 


SUBROUTINE    EP.EZE ( PL  0 AC, NO  INT, NY  I  NT ,FMAX ,PMIN, VMAX ,VMAXX) 
SUBROUTINE    CALCULATES    THE    MAXIMUM   AND    MINIMUM    PCkER    PECUIREHENTS 
AND    MAXIMUM    KIND    SPEED 

IMPLICIT    REAL*8IA-H,0-Z> 

REAL'S    PLCAD(8,4>,VMAX1S,41 

PMAX=0.000 

VMAXX=0.000 

PMIN=1.0013 

DO    10   J=1,NVINT 

DO    10    I=1,NDINT 

A=PLCADI I,JI 

B=VKAX(I,JI 

IFtA.GT.PMAXI    PMAX-A 

IF  t  B    .GT.    VMAXXI    VMAXX-8 
0    IFIA.LT.PMINI    PNIN»A 

RETURN 

END 


0001  SUBROUTINE    Mi'NEY  I  PRATED,  VRATECVREF  ,      ANCCST, CM.CRF , ICOST) 

C*»*    SUBROUTINE    CALCULATES   COST    Op    NTGS    (CAN    ACC    AOCITICNAL   COST    FUNCTION    IF 
C  DES1REC 

0002  IMPLICIT   REAL*8(A-H,0-ZI 

0003  DIMENSION    VREFI2I 

0004  REAL'S    C0EF<4,2l/22S7.8CO,-0.465782C0,7.7  39710O,2.57327O-2, 
1    4*0. ODO/ 

0005  IF    (PRATED. GE. 1.0001    GO    TO    10 

0006  UNTCST«CCEFU,ICCST>»PRATE0>»C0EF(2,ICOST) 

0007  GO    TO    20 

0008  10    PRATLN=DLOG(PRATEDI 

0009  UNTCST.DEXP(CCEF(3,ICCST>»CCEF(2,IC0STI»PRATLN»C0EF(4,IC0STj« 
X      PRATLN**2) 

0010  20    CCST»FRATEC«UNTCST»IVREFI ICCST1/VRATECI  »»2 

0011  ANC0ST«C0ST«CRF 

0012  CM=0.03DO«CCST 

0013  RETURN 

0014  END 
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C 

0001  FUNCTION    C-LCUAC    U,B,FN) 

C*«*  GfiUSS-LEGE\DRE    CUADRATURE    CF    FUNCTION    FN    CVER    INTERVAL    1  A, B) 

C***  INTEGRAL     IS    SET    TC    ZERO    IF    LOWER    LIMIT    LARGER    THAN   UPPER    LIMIT 

0002  IMPLICIT    R?AL*3(A-H,C-Z1 

0003  REAL*8  ROOT (201 . WEIGHT(20I 

0  004  COMMON/ L 1NK3/RC0T, HE IGHT»NHALF 

0005  GLQUAC=O.OD0 

0006  IF  (A.GE.51  RETURN 
000?  BA=0.5D0*(B-A1 
00C8  AB  =  C.5DC»(A*-B) 

0009  DO  10  I=1,NHALF 

0010  10  GLQUAD=GLJUAD*-NEIGHT(  I J  *l  FN(  AE*BA*RQOT(  I  l)  +  FN(  AB-SA*RGOT(II  11 
0  011  GLCUAD=eA*GLgUAO 

0012  RETURN 

0013  ENO 


0001  FUNCTION    FHV1 

C***  SUBROUTINE  COMPUTES  VALUES  CF  EITHER  { V/VRATEO ) **3»BETA  OR  JUST  BETA 

C  DISTRIBUTION 

0002  IMPLICIT    PEAL*e(A-H,G-Z) 

0003  C0NHCN/L!NK2/AU,BIJ»VHAXIJ,VRATED,  PRATED, PRACT3»FRACT4 

0004  FI=(V/VMAXIJ)**IAIJ-1.DO)*(1.CO    -V/VMAXI J J**(BJ J-l.CO    ) 

0005  RETURN 

0006  ENTRY  V3FI{V1 

0007  V3FI=(V/VRATEDJ**3*(V/VMAXIJ]**IAIJ-1.D0)*(1.D0-V/VHAXIJJ**IBIJ- 
1  1.001 

0008  RETURN 

0009  END 
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00000300 

TO    FIND    THE    UNCONSTRAINED    MINIMUM    CF    A    FUNCTION   OF    MANY 

00300500 

VARIABLES    BY    SIMPLEX   PATTERN    SEARCH   METHOD    STARTING   FRCM 

DOD00600 

AN    ARBITRARY    PCINT    ENTERED. 

OOOCO700 

oooooaoo 

00000900 

DESCRIPTION    OF    PARAMETERS 

00001000 

FX             -    WITH    (Nl    DIMENSIONS,    THE    ENTERING    PCINT    WHEN 

00001100 

CALLING    AND   THE    MINIMUM    PCINT    WHEN    RETLRNING. 

00001200 

FY              -    FUNCTION    VALUE    AT    RETURNING    PCINT    IN    RETURN. 

00001300 

N                -    NUMBER   OF    VARIABLES    CF    THE    PROBLEM. 

00001400 

0                -    WITH    (Nl    DIMENSIONS,    THE    STEF-SIZES    FOR    EACH 

00001500 

OIHENSICN    FOR    INITIAL    SIMFLEX    SET-UP. 

00301600 

ITOUT      -    FREOUENCY    CF    W1THI N-SEARCH    INTERMEDIATE    PRINT-OUT 

00001700 

DESIRED.    WHEN    PUT    1T0UT=THE    NUMEER    IN      ITMAX, 

00001800 

NO    INTERMEDIATE    PRINT-OUT    WILL    EE    OUTPUT. 

00001900 

ITMAX      -    MAXIMUM    NUMBER   SEARCH    ITERATION    ASSIGNEO.    WHEN 

OOOC2000 

EXCEEDED,    THE    SEARCH    WILL    EE    TERMINATED    AND    RETURN 

03002100 

THE    LAST    MINIMUM    DATA    SEARCHED. 

00002200 

EPSI         -    STOPPING    CRITERION.    WHEN    EPSI    -CE.    SY,    STANDARD 

00302300 

DEVIATION    CF    FUNCTION    VALUES    EVALUATED   AT    CURRENT 

00002400 

SIMPLEX    VERT  ICES, RETURN    THE    MINIMUM    DATA    . 

00002500 

ALPHA      -    REFLECTION    COEFFICIENT,    SUGGESTED    VALUE    IS    1.0    . 

00002600 

BETA         -    CONTRACTION    COEFFICIENT,     SUGGESTED    VALUE    IS    0.5    . 

00002700 

GAMMA      -    EXPANSION    COEFFICIENT,    SUGCESTEC    VALUE    IS    2.0    . 

000C2800 
00002900 

SUBROUTINE    NEEDED 

OC003000 

SUBROUTINE    0BJN(X,Y,N1         -    FOR    COMPUTE    FUNCTION   VALUE    Y 

00003100 

*T    Y!I  »  ,    W^EPE;           1  =  1.?:           :N     . 

OO003300 
00003400 

OOOl 

SUBROUTINE    S IMPXIFX ,FY,N , D, I  TOUT, ITMAX.EPSI ,ALFHA,BETA, GAMMA! 

00000400 

000035CO 

0002 

IMPLICIT    REAL*8(A-H,0-Z) 

) 

0003 

REAL*E    XI9,8),YI9I ,FX(N».0iN) 

00003600 

0004 

1003 

FCRMAT(3X,fc5(lH*M 

00003700 

0005 

1011 

FORMATI 5X.5H0Y    =    ,E11.5,9H      ITER    -    ,!4,1CH      NOFT.    »    .I4.10H      NOCVN00003800 

I    «    ,141 

00CC39O0 

0006 

1012 

F0RMAT(7X,8HN0RFT    -    ,14 ,4X, 8HN0EXP    •    ,14,10H      NCCNT    =    ,14,10H      NOC00004000 

1UT    «    ,141 

00004100 

0007 

1C13 

FCRMATI7X.24HCURRENT    SEARCHED    DATA    ../10X.3FY*    ,E11.5,1H.I 

00004200 

9008 

1014 

FORMAT(10X,2HX{  ,13, 4H)    =     ,£11  .5  ,1H  ,  .5X.3HCX  ( ,  13  ,4H)    =    ,EU.5,1H,1 

00004300 

0009 

1015 

FCRMAT(7X,8HYMEAN    =     ,E15.E,9H    ,       SY    *    ,El!.8,2h    .i 

00004400 

0010 

1016 

F0RHAT(5Xi24H**CtJT    STEP-SIZES    TIMES       ,I3,2H    .) 

00304500 

0011 

1023 

F0RMAT(5X,26H»»    ITERATION   NO.    EXCEEDED    ,I5,2H    .) 

00004600 

0012 

KULT-1 

O0C047C0 

0013 

NOPT«0 

00004800 

0014 

NDCUT«0 

00004900 

0015 

NCCVN-0 

00005000 

0016 

I1ER-0 

00005100 

0017 

NORFT-0 

00005200 

001! 

NOEXP-0 

00005300 

0019 

NOCNT-6 

00005400 

0020 

FN=N 

OC0C55O0 

0021 

NN=-N*1 

00005600 

0022 

CALL    CBJN(FX,YF,NI 

00005700 

0023 

LCCAT*! 

C0005800 

002* 

IWAY»1 

00305900 
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C         SET    UP    INITIAL    SIMPLEX    . 

00006000 

0025 

2 

00   6    J=1,N 

00006100  ■ 

0026 

00    3    1-l.J 

00006200 

002T 

3 

X(I,JI=      FX{J»-D(J> 

00006300 

002B 

FJ»J 

00006400 

0029 

XUMiJI*      FX(J)»FJ«0(JI 

00006500 

0030 

IFI J-N)4,6,6 

00006600 

0031 

4 

JH*J*2 

00006700 

0032 

DO    5    I*JM.NH 

00006800 

0033 

5 

XII, J>-      FXIJI 

00006930 

0034 

6 

CONTINUE 

00007000 

0035 

00    8    1=1. NN 

00007100 

0036 

DO    7    J«1,N 

00007200 

0037 

7 

FX(JI'X(I>JI 

00007300 

0038 

CALL    C6JNCFX,YF,NI 

00007400 

0039 

8 

YIII-YF 

00007500 

0040 

INI  =  1 

00007600 

C         REARRANGE    ORDER    (OVERALLI    . 

00007700 

0041 

9 

I"l 

00007800 

0042 

NS=NH 

00007900 

0043 

10 

IF(Y(I)-r{NS>)13,U.lI 

00008000 

0044 

11 

YTE*«YCNSI 

00008100 

0045 

YINS)=YCI) 

00003200 

0046 

Y(1I»YTE« 

00008300 

0  04T 

DO    12   J-l.N 

00008400 

0048 

FXIJ)-XiNS,JI 

00008500 

0049 

X(NS,  Jt-XU.J! 

00008600 

0050 

12 

XU»JI*      FX(J) 

00008700 

0051 

13 

IFINS-I-1) 15,15,14 

00008800 

"052 

14 

NS=N5-I 

CO0CS5CC 

0053 

GO    TO    10 

00009000 

0054 

15 

I«t*l 

00009100 

0055 

IF(!-N-1I16,17,17 

00009200 

0056 

16 

NS«Nfl 

00009300 

0057 

GO   TC    10 

00009400                    ) 

0058 

17 

IFdNII    65,65,501 

00009500 

0059 

501 

LGCAT=2 

00009600 

0060 

IWAY=2 

00009700 

0061 

GO    TO    120 

00009800 

C         CCKPt     ■ 

00009900 

0062 

16 

DC    20   J=1,N 

00010000 

0063 

PXT»X(1,JI 

00010100 

0064 

00    19    1=2, N 

00010200 

0065 

19 

PXT=PXT+X(I ,J) 

00010300 

0066 

20 

X(W-2,JI  =  PXT/FN 

00010400    ■ 

C         »• 

."AKE    REFLECTION    HOVE    . 

00010500 

0067 

DO    21    J«1,N 

00010600 

0068 

X(N+3,JI  =  XIN+2,J»  +  ALPHAMXIN«-2,JI-X(N+1,J)I 

00010700 

0  069 

21 

FX( JI=XtN+3,JI 

00010800 

0070 

CALL    CSJN(FX,YF,N) 

00010900 

0071 

YCN»3I=YF 

00011000 

0072 

N0PT=N0PT+1 

C0011100 

0073 

LCCAT»3 

00011200 

0074 

IWAY*3 

00011300 

00'5 

GO   TO    500 

00011400 

0076 

22 

1F(Y(N«3I-YI1II29,23,23 

00011500 

0077 

23 

IF(YIN*3)-Y(N>]24,26,26 

00011600 

0  078 

244 

IWAY«7 

00011700 
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0079 

24    DO   25    I»1,N 

00011800 

0080 

25    X(N*1,I)=X(N+3,II 

00011900 

0081 

Y(N*1I»Y(N»3I 

00012000 

0082 

ITER=ITER»1    . 

00012100 

0083 

NCRFT«NORPT*l 

00012200 

0064 

GO   TO    100 

0C012300 

0085 

26    IFtY(N+3>-Y[N*ll>27,49,49 

00012400 

0086 

27    DO    28    1-1. N 

00012500 

0037 

28    X(N<-1,I>=X(N»3,I) 

00012600 

0088 

Y(NH>=Y(»»3) 

00012700 

0  089 

I7ER=ITER*1 

00012800 

0090 

NDRFT*NORFT*l 

00012900 

0091 

GO    TO  49 

00013000 

C         **HAKE    EXPANSION    HOVE    f 

00013100 

0092 

29    DC  30    J  =  1,N 

00013200 

0093 

X[NM.JI  =  X(N+2,JKGAKMA*(X(N+3,JI  -X[N*2,J1) 

00013300 

0  09* 

30    FXIJI«XIN»4,J> 

00013400 

0C95 

CALL    CBJMFX,YF,N1 

00013500 

0096 

Y(H*4)*YF 

00013600 

0097 

NCPT-N0PT*1 

00013700 

0098 

LOCAT-4 

00013800 

0  099 

IWAY-4 

00013900 

0100 

GC   TO    500 

00014000 

0101 

31    IF(Y[N<-41-Y(1>>32,244,244 

00014100 

0102 

32    DO    33    1=1, N 

00014200 

0103 

33    X<N»1,1)«XI(H-*,I1 

00014300 

0104 

Y(N*11-Y1N»4I 

00014400 

0105 

ITER=ITERH 

OOC14500 

0106 

NCEXP=NOEXP*i 

00014600 

0107 

GO   TO    100 

onol47nn 

C         "MAKE    CONTRACTION    MOVE    . 

00014300 

0108 

49    DO    50    J-l.N 

00014900 

0109 

XfN+5,J»=X(N+2,J]+BETA*{XtN+l,J)-XIN+2,J)l 

00015000 

0110 

50    FXI JI=X(N+5.JI 

00015100 

0111 

CALL    C8JN(FX,YF,N) 

00015200 

0112 

Y(N+5)=YF 

00015300 

0113 

N0PT*N0PT+1 

00015400 

0114 

LOCAT-5 

00015500 

0115 

IhAY-5 

00015600 

0116 

GO   TO    500 

00015700 

0117 

51    IFIY<Nt51-YINH)>52,6C,60 

00015800 

0118 

52    DC    53    I-l.N 

00015900 

0119 

53    X(N*1,II=X(N*5,1I 

0D016000 

0120 

YIN»1).Y(N«5I 

00016100 

0121 

ITER=ITER»1 

00016200 

0122 

N0CNT=NCCNT*1 

00016300 

0123 

NCCVN=N0CVN*1 

00016400 

0124 

GO   TO    110 

00016500 

C         «*CUT    DOWN    STEP-SUES    . 

00016600 

0125 

60   00    62    1=2, NN 

00016700 

0126 

DO    61    J=1,N 

00016800 

0127 

XII,JI<[Xll,JI«X<I,JJ)/2.OD0 

00016900 

0128 

61    FX(J)=X(I,J1 

00017000 

0129 

CALL    CBJN(FX,YF,N) 

00017100 

0130 

62    Y( I)*YF 

00017200 

C         "REARRANGE    ORDER    (OVERALL    )    . 

00017300 

0131 

IN!=0 

00017400 

0132 

GO    TO   9 

00017500 
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ABSTRACT 

To  obtain  accurate  estimates  of  extractable  energy  from  the  wind  at 
a  given  location,  it  is  first  necessary  to  obtain  an  accurate  description 
of  the  expected  distribution  of  wind  speeds.   For  such  analyses,  the  use 
of  analytical  representations  for  the  wind  speed  distributions  often 
simplifies  the  calculation  as  well  as  smooths  out  statistical  fluctuations 
in  the  observed  wind  speed  data.   In  the  first  phase  of  this  study,  three 
techniques  for  fitting  a  Weibull  distribution  to  observed  wind  speed  data 
are  examined.   In  addition,  the  beta  distribution  is  introduced  as  an 
alternative  wind  speed  distribution  model  and  a  mat chin g-moments  scheme 
is  presented  to  obtain  the  beta  distribution's  parameters.   Two  goodness 
of  fit  tests  are  performed  on  each  analytical  distribution  to  test  the 
appropriateness  of  each  model  in  describing  544  observed  wind  speed 
distributions.   It  was  found  that  least  squares  fitting  techniques  pro- 
duce Weibull  distributions  which  poorly  represent  wind  data,  but  that 
both  Weibull  and  beta  distributions  give  excellent  fits  to  the  data  when 
the  parameters  are  obtained  with  a  matching-moments  technique. 

In  the  second  phase  of  this  work,  the  analytically  fit  wind  speed 
distributions  are  used  in  a  methodology  to  select  the  optimally  sized  wind 
turbine  generator  system  for  given  demand  power  requirements.   Such  an 
optimally  sized  system  will  yield  the  maximum  net  economic  savings  for 
an  enterprise  which  uses  the  system.   In  particular,  the  sensitivity  of 
the  optimal  wind  system  to  various  problem  parameters  such  as  electricity 
cost  and  wind  and  load  characteristics  are  investigated.   It  was  found 


that  to  accurately  estimate  the  capability  of  a  particular  wind  turbine 
system  to  supply  the  energy  needs  of  a  specified  demand  load,  daily  and 
seasonal  wind  speeds  and  demand  loads  as  well  as  diurnal  variations  in 
these  characteristics  must  be  known. 


